Unbalanced Optimal Transport through Non-negative Penalized Linear Regression
Authors: Laetitia Chapel
Rémi Flamary
Haoran Wu
Cédric Févotte
Gilles Gasso
Venue: NeurIPS 2021 Poster (2021)
Keywords: Optimal Transport
Unbalanced Optimal Transport
Penalized Linear Regression
Lasso
MM algorithms

Abstract: This paper addresses the problem of Unbalanced Optimal Transport (UOT) in which the marginal conditions are relaxed (using weighted penalties in lieu of equality) and no additional regularization is enforced on the OT plan. In this context, we show that the corresponding optimization problem can be reformulated as a non-negative penalized linear regression problem. This reformulation allows us to propose novel algorithms inspired from inverse problems and nonnegative matrix factorization. In particular, we consider majorization-minimization which leads in our setting to efficient multiplicative updates for a variety of penalties. Furthermore, we derive for the first time an efficient algorithm to compute the regularization path of UOT with quadratic penalties. The proposed algorithm provides a continuity of piece-wise linear OT plans converging to the solution of balanced OT (corresponding to infinite penalty weights). We perform several numerical experiments on simulated and real data illustrating the new algorithms, and provide a detailed discussion about more sophisticated optimization tools that can further be used to solve OT problems thanks to our reformulation. 

Unbalanced Optimal Transport through Non-negative Penalized Linear Regression
Laetitia Chapel'
Rémi Flamary" CMAP, Ecole Polytechnique
IRISA, Université Bretagne-Ste SSd
Vannes, France laetitia. chapel@irisa
Palaiseau, France
remi lamaar@@polytechniuueeeeeo edu
Haoran Wu LITIS & IRISA Rouen & Vannes, France haoran. vu@univ-ubs frr
Cédric Févotte
IRIT, Université de Toulouse, CNRS
Toulouse, France
feddii..fvvotteirit..
Gilles Gasso
LITIS, INSA Rouen Normandie Rouen, France
Abstract
This paper addresses the problem of Unbalanced Optimal Transport (UOT) in which the marginal conditions are relaxed (using weighted penalties in lieu equality) and additional regularization is enforced the OT plan. In this context, we show that the corresponding optimization problem can be reformulated non-negative penalized linear regression problem. This reformulation allows us to propose novel algorithms inspired from inverse problems and nonnegative matrix factorization. In particular, we consider majorization-miiiiizztoooo which leads in our setting to efficient multiplicative updates for variety of penalties. Furthermore, we derive for the first time an efficient algorithm to compute the regularization path of UOT with quadratic penalties. The proposed algorithm provides a continuity of piece-wise linear OT plans converging to the solution of balanced OT (corresponding to infinite penalty weights). We perform several numerical experiments on simulated and real data illustrating the new algorithms, and provide detailed discussion about more sophisticated optimization tools that can further be used to solve OT problems thanks to our reformulation.
Introduction
Optimal Transport (OT) theory provides powerful tools for comparing probability distributions and has been successfully employed in wide range of machine learning applications such supervised learning (Frogner et al. 2015), clustering (Ho et al., (017)). generative modelling (Arjovsky et al. 2017). domain adaptation (Courty et al. 2017). learning of structured data (Maretic et al. 2019 Vayer et al. 2019) or natural language processing (Kusner et al.) 2015). among many others. One reason for those recent successes is the introduction of entropy regularire OT that can solved with the efficient Sinkhorn- Knopp matrix scaling algorithm (Cuturi 2013). However, the classical OT problem seeks the optimal cost to transport all the mass from source distribution target one
First two authors have equal contribution
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

Viillani 2009), , greatly limiting its use in scenarii where the measures have different masses or when they contain noisy observations or outliers.
Unbalanced Optimal Transport (UOT) (Benamou 2003) has been introduced to tackle this shortcom- ing, allowing some mass variation in the transportation problem. It is expressed as a relaxation of the Kantorovich formulation (Kantorovich] 1942) by penalizing the divergence between the marginals of the transportation plan and the given distributions. Several divergences can be considered, such as the Kullback-Leiber (KL) divergence (Frogner et 2015 Liero et al. 2018)) the & norm corresponding to the partial optimal transport problem (Caffarelli and McCann 2010 Figalli 2010). or the squared lg norm (Benamou 2003). Regarding numerical solutions, CChzzt et al. 2018) considered an entropic- regularized version of UOT leading to a class of scaling algorithms the vein of the Sinkhorn -noppp approach (Sinkhorn and Knopp [967). The introduction of this entropic regularization improves the scalability of OT, but involves spreading of the mass and a loss of sparsity in the OT plan. When a sparse transport plan is sought, the convergence is slowed down, necessitating the use of acceleration strategies (Thibault et al. 202T). Regarding UOT with the (squared) Eg norm, Blondel et al. (2018) showed that the resulting OT plan is sparse and proposed to use an efficient L-BFGS-B algorithm (Byrd et al. (995) to address this case. Note that the L-BFGS-B method can be used to solve UOT with differentiable divergences even without the entropic-reegulrizatioo the Sinkhorn-like iterations. Also note that, as for balanced OT, UOT can be solved more efficiently when the data has specific structure, such unidimensional distributions (Bonneel and Cocurjolly) 2019) or distributions supported on trees (Sato et al. 2020). Finally, recent work investigated UOT between Gaussians and provided closed form solutions for the regularized Janati et al. (2020) and unregularized (Janati 2021 Eq. 2.72) versions of UOT associated with a KI. divergence.
the OT plan that induces
Contributions. In this paper, we show after some preliminaries that UOT can be recast as a convex penalized linear regression problem with non-negativity constraints (Section 2.2) The main interest of this reformulation resides in the fact that non-negative linear regression has been extensively studied in inverse problems and machine learning, offering a large panel of tools for devising new numerical algorithms. Our reformulation involves a design/dictionary matrix that is structured and sparse. Leveraging this structure, we propose two new families of algorithms for solving the exact (i.e., without regularization of the plan) UOT problem in Section 3
We first derive in Section 3.1 a new iijirizattio-Mininiizz with Bregman divergences, and more specifically KL and 6z-penalized UOT. The MM approach results in multiplicative updates that have appealing features: i) they are easy to implement, ii) have low complexity per iteration and can be instantiated on GPU, iii) ensure monotonicity of the objective function and inherit existing convergence results. Our methodology is inspired by well-known algorithms in image restoration (Richardson) 1972 De Pierro| 1993) and non-negative matrix factorization (NMF) (Lee and Seung 2001] Dhillon and Sra [2055 Févotte and Idier| 2011). Interestingly, the resulting multiplicative updates bear a similarity with the celebrated Sinkhorn scaling algorithm, with some key differences that are discussed.
(MM) algorithm for solving UOT
Next, we derive in Section 3.2 an efficient algorithm to compute the regularization path in é2- penalized UOT. To do so, we build on our proposed reformulation and more precisely on the fact that fz-penalized UOT can be reformulated as a weighted Lasso problem. We propose new methodology inspired by LARS (Efron et al. [004] Hastic et al., 2004). which, to the best of our knowledge, is the first regularization path algorithm for OT problems. It brings novel understanding of the properties of the evolution of the support of OT plans, besides the practical interest of computing the complete regularization path when hyperparameter validation is necessary.
Our new families of algorithms (MM for general UOT, LARS for é2-penalized UOT) are showcased in the numerical experiments of Section 4 Python implementation of the algorithms, provided in supplementary, will be released with MIT license on GitHub. The connection between UOT and linear regression that we reveal in the paper opens the door further fruitful developments and in particular to more efficient algorithms, thanks to the large literature dealing with gon-negative penalized linear regression. We discuss those possible research directions in Section 3 before concluding the paper.
Notations. Vectors such are written with lower case and bold font, with coefficients m; or [m]; according to context. The [A)-dimensional sub-vector with indexes in set A is written m,A Matrices such as M are written with upper case and bold font, with coefficients Mij. We introduce
vectorization operator defined by m vec(M) [M M1.2.... ... Mp 1,, Mn. i.c., the

---

concatenation of the rows of the matrix, following the Numpy/C memory convention. 1n is a vector of n ones and M 2 0 denotes entry-wise non-negativity. Finally, D the Bregman divergence generated by the strictly convex and differentiable function y, i.e., (u, =,, (ds, = [x()((u) - (p() (w; [=]]-
Reformulation of UOT as non-negative penalized linear regression
2.1 Background on Optimal Transport
Let us consider two clouds of points X = (x= and Y = (y,) Let E R and b E IRj be two discrete distributions of mass on X and Y. such that (resp. 6,) is the mass at æ; (resp. y,). The balanced OT problem, as defined by Kantorovich (1942), is linear problem that computes the minimum cost of moving to b:
(T b) min CT such that (s.t.) T1m a, T,T 1n =b
(1)
where (.-) is the Frobenius inner product, T E Rxxm the transport plan and C B the cost matrix. The entry Ci,j of C represents the cost of moving point y)- The Wasserstein 1- distance (also known as the earth mover's 's distance) is obtained for Ct,- = || Výll- The constraints on the transport plan T require that ||a||] ||b]]] and that all the mass from is transported to b. These constraints can be alleviated through relaxation, leading to UOT (Benamou 2003):
UOT (a, b)
min T20
(C,T) + DD a) + >2D In, b).
The deviations from the true marginals are penalized by means of a given Bregman divergence Dip. as introduced in Chizat et al. (2018), where a₁ and X2 are hyperparameters that represent the strengths of penalization. Note that, in the case of ||a||]
[[b]] balanced OT (Eq. is recovered when A1 = X2 + oc. Furthermore, when a or x2 > co, we recover semi-relaxed OT (Rabin et al. 2014). In practice, authors often set d A2 > for UOT in order reduce the necessity of hyperparameter tuning. Various divergences have been considered in the literature. The ε norm gives rise so-called partial optimal transport (Caffarelli and McCann 2010). The squared norm provides a sparse and smooth transport plan (Blondel et al. (2018) when introducing a strongly convex term in Eq. (2) Chizat et al. (2018) derive efficient algorithms solve Eq. (2) for several divergences by adding an additional regularization term >reg (T (T, ab ). In particular, entropic regularization is obtained when the KL divergence is used, promoting a dense transport plan unlike exact UOT.
2.2 Reformulation of UOT
UOT cast as regression. Let t = vec(T)) c =ve(CC) and T (a. b"). Problem 2 can be re-written as
(t) def X D.
min 220
(3)
and as such be expressed as non-negative penalized linear regression problem, where the design matrix H = H,, H. IT is the concatenation of the matrices H, and H. that compute sums of the rows and columns of T. respectively (see expressions in Section ?? of the supplementary material). Note that, for the sake of simplicity, we consider here A X2 = À but this hypothesis could be easily alleviated for a given family of divergences (see Sec. 5 for a discussion). Important features of Eq. 3 should be discussed. First, (t) is convex thanks to the convexity of Bregman divergences w.r.t their first argument. Second, H is very structured and sparse (with a ratio of only m+n non-zero coefficients) which will allow for more efficient computations and updates than with dense H. Finally, since & 2 0 and 2 0, the linear term can be expressed as _ t= E. cat. = 1 Σ. Cslal- This corresponds
weighted & regularization, promoting sparsíty in and hence in the transport plans. Note that the "sparse" regularization is here controlled by & (instead of À in classical penalized linear regression), meaning that the sparsity promoting term will be more aggressive for small A.
Solving problem 3 Problems of the form of Eq. are well-known in inverse problems and NMF In inverse problems, typically acts as a clean image degraded by operator H (e.g., a convolution) and noise. The data fitting term D. (tt, z) captures assumptions about the noise corrupting the observed
3

---

image y. Sparsity is common regularizer of t. In NMF, given set of nonnegative samples (yr one wants to learn a non non-negative dictionary H and non-negative owver-dimessionll embeddings {tr) such that Y₁ R H4 (Lee and Seung 1999). Updating the latter involves optimization problems of form (3) (with or without sparse regularization).
contrast to problem (3). the data fitting term more commonly (y, Ht. instead of D. (Ht, y) in inverse problems and NMF. This is because the former is a log-lkkelihood in disguise for the mean-parametrizee exponential family, and takes important noise models as special cases, such as Poisson, additive Gaussian or multiplicative Gamma noise Tééotte and Idier (01T)) Using such penalizations with reversed arguments would be possible in our 2018) for simplicity.
well but we stick to the now standard formulation of (Liero al., 2018 Chizat et al.
In the next section, we will first leverage classical family of algorithms in inverse problems and NMF, namely MM, to obtain new algorithms for KI. and fz-penalized UOT (possibly with entropic regularization in the first case). Second, we will leverage results about non-negative Lasso to design an efficient algorithm to compute the regularization path of fz-penalized UOT.
Novel numerical solvers for UOT
3.1 Majorization-MMinnzzztii (MM) for UOT
General MM framework. MM algorithms have been around long time in inverse problems and
NMF to solve problems of form Classical algorithms for NMF such (Lee and Seung. 2001)
have built on seminal MM algorithms for inverse problems such as (Richardson) 1972. De Pierro
[993). Subsequent works NMF such as (Dhillon and Sra 2005 Févotte and Idier 201 Yang and
201T) have further contributed novel MM algorithms for larger classes of problems, including
larger families of divergences. In nutshell, MM consists in iteratively building and minimizing an upper bound of the objective function which is tight at the current parameter estimate (and referred to as auxiliary function), see Hunter and Lange (2004); Sun al. (2017) for tutorials. In NMF, common approach consists of alternating the updates of the dictionary H and of the embeddings. In our case, H is fixed and we may use the results of (Dhillon and Sra. 2005) build an auxiliary function for term D& Ht, y), to which we may simply add the linear term c t// to obtain valid
auxiliary function for (t). Let i denote the current estimate of t, 7 Σ, H.
H., and
Gx(l,t) En.0
Higst Zii
Gj
tú cst,
(4)
2
where cst EEle'(n) - sp(vl) Then, Gx t) is an auxiliary function for Ex(t), i.e. Vt. &x(t,) > Ft(t) and G2 2) Fx)... Let (k-1) G2 t(K)), then Fx (t(k)) () ((0)) 2 Gx(e(++1), (())) 2 (x((()))), producing descent algorithm over F. The trick to 2.,
Hiit,), 2.,, thanks to the convexity of yo, see details in (Dhillon and Sra 2005). We provide below the resulting algorithms for the KI. and lg penalizations, with detailed computations available in Section ?? of the supplementary.
obtain Gs to apply Jensen inequality to ((E, Higt)) <tE, Zr.5 Meit,) S E3 Znjip(
MM for KL-penalized UOT. The KL divergence is obtained with sp()) = y log y y. Minimizing G (t, (()) in that case leads following multiplicative update:
[
log(y) .H log [Ht (b) I, to)
((++1)
(00)
cxp
(5)
H 11
Owing to the structure of and H, , the update can be re-written in the following matrix form:
( diag
(T())) o exp
C
diag
b
(6)
(
2x
T
where @ entrywise multiplication and divisions are taken entrywise as well. The multiplicative
update 6 is remarkably similar to the well-known Sinkhorn-Knop algorithm that has been used in

---

numerous OT problems involving KL regularization. But instead of two separate steps for the left
and right scaling, Eq. (6) applies these scalings simultaneously in a unique update using the diagonal matrices (and a form of geometrical average). Also note how the scaling factor exp along iterations the coefficients of the transport plan with large costs.
C penalizes 2x
MM for lz-penalized UOT. The quadratic loss is obtained with plv) = 2/- In that case, minimiz- ing Gx(t, (()) s.t. non-negativity leads following multiplicative update:
max (0, al m 1 2C)) T O T,T
((++1) (h)
with 0z= L
(7)
Interestingly enough, update 7 prunes any coefficient Tij in T such that a; bb - C < 0 from the very first iteration, providing a useful certificate on the support of the solution.
3.2 Regularization path for Ey-penalized UOT
Let us focus on the case where D, is quadratic divergence. As mentioned in Section 2.2 Eq. is then positive weighted Lasso problem, allowing us to derive the first regularization path algorithm for computing the whole set of solutions for a varying À from 0 to +o0. Note that the path 's extreme point recovers the balanced OT solution. We show that the path is piecewise linear 1/A between changes in the active set
supp(t*), where &* vec(T^) and T is the OT plan for given hyperparameter A. The main steps of the algorithm are roughly as follows: given current solution (Ag, Tx and current active set Ag, look for the next value dR+1 > As such that the active set changes (i.c., A&+1 * Ax), either because one component enters or leaves the active set. We describe our algorithm below.
KKT conditions of the Ez-penalized UOT problem. The Lagrangian for problem B
writes:
1 1 c + Ht (H -y) - (t, =
(8)
X 2
where " represents the Lagrange parameters. We denote m = H y = vec(al, + 1,bT). KKT optimality conditions state that i) Ve,x A) = c H' Ht m -y (stationarity condition), ii) y G t (complementary condition) and iii) Y 2 (feasibility condition).
Piecewise linearity of the path. Assume that, at iteration k, we know the current active set A = A and we look for ti (the other values of L being 0). Let HA. mA and Ca denote the corresponding sub-matrix and vectors (see Appendix ?? for rigorous definitions). Because of the complementary condition, we have Y O. Using À A + e. with e > 0 small enough ensure that the active set remains the same, the stationarity condition writes:
HHH = mA 7 CA
+
ta mA
Eq. 9 shows that the optimal ta (and hence tー) can be solved for any X E [Ak, A6411, i.e., when the active set A remains the same, by solving linear problem. It also reveals the piecewise linearity in A-1 of the path when is fixed. As expected, balanced OT is recovered when À 00.
with m = (HAHHHHH)--m
and ēA (HH
Finding A given As). Given a current solution (Ae, t>a) and X d c, we increase the until we reach a change in the set of active components. This happens whenever the first of the following two situations occurs.
One component in A becomes inactive. In that case, we remove the index E A with the smallest X Dy. that violates the constraint. In such case, [m = [č,]i/A and may write
ČA
min
(10)
>>
m
where ming A, indicates the minimum value in the vector greater than d and the division is entrywise.
5

---

Algorithm I Regularization path of ee-penalized UOT
Require: a, b, C, A 0, to 0, A
Ao
0, =1
A1 = min CA A = A = arg min
CA.
H H =2
m
m
Ai md
44
A1
while ( - y) (Ht y) + 0 do
Xy. da f Compute as in Eq. (10) and Eq.
AR+1 - min(A, )a)
At+1 H - - H
A
MA
A Ak+1 f- Update active set for next iteration.
H A) + Update from (H). HA,) with Schur complement (see supplementary ??)
1 + +1
end while
return (Ak, )x
One component in A becomes active. This occurs when the KKT positivity constraint V becomes violated. Assume this happens at index i E A for the smallest value Xa > d of X. In such case, the stationarity condition outside the active set be rewritten:
1 }0)] ma Hea+ H ( +
CA
H
HC H
-Ivali +
Aa min >>.
(11)
mA
H
where mn (resp. č) equals ma (resp. ε사) on A and zero on A.
In practice, at each step of the path, compute both X, and X and set Ak+1 min(Ar, Aa). The active set As+1 is obtained by either removing the index i E A correponding of the arg min of eq.
(10) (case AA+11 X,) by adding the index E A corresponding to the arg min of eq. (II (case Ak+1
A,).
involves the computation of the matrix (HAHA) which is of size |A| x [AA. As only one index leaves or enters the active set at each iteration, we can use the Schur complement of the matrix to compute its value from (H X, H AA) alleviating the computational burden of the algorithm (see Section ?? of supplementary). Algorithm I sums up the different steps of the full path computa- tion. At each iteration, we compute Au, 'y, update the inverse matrix (H A₁ H A.) and estimate the solution A with complexity of O(m)).
Numerical computation of the entire path. Eq. 9
only involves matrix-vector computations
Regularization path of the semi-relaxed & le-penalized UOT. As a side result, let consider the semi-relaxed OT problem SROT* (a. b) = minT22,T =b (C,T) + Al|TL,.. m a The main difference with UOT that the equality constraint T I, b (equivalent to H.t b) must always be met. This leads to the following Lagrangian:
Lx 7, = Tt t [ (Ht (H,-aj) a) (H, - a) + (H,t b) u - Tt,
(12)
where u E I contains the Lagrange parameters associated to the equality constraints. The KKT optimality conditions now dictate that i) Lz u) - H H,t H.a a H.- ug-u-a O. ii) y 3 0, iii) 7 2 and H b O. We can use the same reasoning than previously to compute the entire path. Details are provided in Section ?? of the supplementary. The main difference lies in solving, at each iteration, a linear system of size (m + Al) to comply with the marginal equality constraint. The path initialized as follows: the jth column of TO for X = 0 is set to the weighted canonical vector by- b es,, where i' argmminCa...b

---

Evolution plan
Cost matrix
OT
0.5
T
0.4
0,3
Gi
T
T22
Th 42
T12
0.2
Tea
T5.1
Tap
7ks
ny 00-03
T30
Tya
6.602 62
0.5
by
A1
Ag
Ag
As
Ag
A7
Figure 1: (Left) cost matrix C' (the higher the cost, the darker the color); (middle) OT plan whose cells are color-coded with respect to the À values which they are activated. The blank cells never enter the active set as the corresponding cost it too high; (right) evolution of Ti when increases. Note that the a-axis is in log scale and discontinued (but still monotonic) between À7 and co.
Distributions
UOT
-
fg
Source
Semi-relssed iy UOT
Semicelaced by A-100
Semiine aced & UDT 1000
Semirelased by UCTT Amind
OT plan
Target marginal Source nangrol
Figure 2: Regularization paths for 2D empirical distributions for Ez-penalized UOT (top) and semi- relaxed UOT (bottom). The OT plan is shown green lines between the source and target samples when Tag > 0 and the resulting marginals are shown as filled circles.
4 Numerical experiments
In this section, we first show the solutions obtained with our solven on simple and interpretable examples. We then evaluate the computational complexity of the different algorithms and finally we show an application where the regularization path can be used on a domain adaptation problem.
Illustration of the algorithms. We first illustrate the regularization path for l (y-penalized UOT on simple example between two distributions containing 3 points each, with different masses and cost matrix C given Fig. [ (left). We can see on Fig. [ (right) that, starting from A = 0 and T= == 0, we successively add remove components in the active set A when increasing the À values. When = oo, we recover the balanced OT solution. Recall that the path is linear in 1/X (and not A). We then illustrate the path for both &2 penalizee UOT and semi-relaxed UOT on two 2D distributions with n = m = 100 samples. We can see in Fig. 2 the difference between the two regularization paths for specific values of A. UOT starts with empty plan for À 0 and then activates samples from both source and target from the closest to the farthest ones until convergence to the balanced OT plan. Semi-relaxed UOT starts with all target samples active due marginal constraints and progressively activates the source samples.
Comparison of the performances of the algorithms. We provide empirical evaluation of the running times of the proposed algorithms, using 2 sets of 10-dimensional points with n = m and drawn according IID Gaussian distributions. The cost matrix C is computed using a squared &2
"Our implementation of the regularization path has been contributed to POT Flamary et al. (021) and the MM algorithms are provided in the repository lattps://github com/lchapell UOT though pettttiee-aiiiedd-eeee Linear- regression
7

---

Regularization path
Eeperalized
KL-penalized
105
10P
10
(nelre
Time (C00)
217
T (twe)
200
10
BrdS
Mulipice updare
107
11000 jelee]
(1500 (CF)
10-1
RPGS
Rague aaccrpat
Mulliposh updre
107
number
Figure 3: (Left) Running times of Alg. 1 w.r.t the number of points: (middle) comparison of ég-penalized UOT with m = 500 (right) likewise for KL-penalized UOT. Dark curves (resp. shaded regions) represent average (resp. variance) values over 5 runs.
norm. We first study the rumning times of the regularization path algorithm, for n n m ranging from 100 to 1000, averaging the results over 5 runs, see Fig. B (left). We empirically observe that log-log plot is near-linear, with an empirical complexity O
in this example.
Using n (Blondel et al. 2018 3 using SciPy (Virtanen et al. 2020) and those of our algorithms: the bz-
500, we compare the running times of the current sfate-of-the-art BFGS algorithm
penalized UOT formulated as a Lasso problem (with both the Celer algorithm (Massias al. 2018)
and the coordinate descent solvers from Scikit-learn (Pedregosa et al. (201)), the multiplicative algorithm for both the &2 and the KL penalties and the regularization path algorithm (see Section ?? of the supplemental material for more details about the solvers and their parameters, together with a comparison of the results of the MM algorithms computed on both CPU and GPU). Figure 3 (middle and right) shows the average running time for all algorithms. For fz-penalized UOT, we observe that, for large A values, the Lasso solvers are the fastest and that, whatever the value X, BFGS is the slowest. We also notice that, for large X, the running times for computing the path remain constant: when the last active set found, computing the OT plan only involves a weighted sum. As for KI.-penalized UOT, the BFGS algorithm is more efficient when large values of À are considered. One can also notice that, similarly to Sinkhorn which is fast for large regularization values, the multiplicative algorithms for both penalties are also fast for high 1/A values.
Regularization path for unbalanced domain adaptation. We demonstrate the interest of having the entire regularization path in classification context where some of the data collection may be polluted by outliers. We consider setup similar to Mukherjee et al. (2020). Let the source X be set of 400 MNIST digits sampled from the digits O. 1, 2, 3 (100 points per class) and let the target Y be set of digits 0, 1 of MNIST (LeCun et al. 2010) and of digits 8, 9 from Fashion MNIST (Xiao et al., 2017). Our setting is simple classification: we classify sample of the target dataset by propagating the label of the source sample provided that the transported mass of the target point is greater than 0.25b,. Note that similarly to Mukherjee et al. (2020) validation set can be used here to select the best A. Figure shows the overall accuracy, defined as the number of samples that are correctly classified divided by the total number of points, and the current accuracy, which is the proportion of well-classified points among the points that classified, i.e., that are receiving mass. One can notice that, as the number of classified points increases (with A), the overall accuracy increases as more and more points are well classified while the current accuracy remains stable until outliers are included in the labeled set. This suggests that UOT can be used not only for classification but also as an automated outlier detection method.
the most transported to,
1.0
0.8
0.6
0.4
0.2
Acc. al Acc. labeled
6.0
100
200
300
100
Numbar laboled points
Figure 4: Evolution of the classi- fication accuracy for the domain adaptation problem w.r.t. the num- ber of classified points.
'Note that we cannot compare our result with solvers for entropic formulation of UOT with Aveg + 0 (that should provide a sparse transport plan) as in that case, their algorithm becomes a fixed point method that cannot converge solution of the problem.
8

---

Discussion and perspectives
We showed that UOT can be recast as a non non-negative penalized linear regression problem, encouraging us dig into this well-established field of research in order to adapt existing algorithmic solutions to the structure of the UOT problem. In this section, we discuss the relation between the proposed algorithms and classical solvers used in OT, and also investigate some research directions that can widen the scope of proposed methods.
Multiplicative algorithms for UOT. As discussed in Section 3.1 the multiplicative updates for the KL divergence obtained from MM resemble the Sinkhorn algorithm from Chizat et al. (2018). except for the joint scaling and the weighting matrix also appears in the Inexact Proximal Point OT (IPOT) algorithm Xie et al. (2020) to solve balanced OT. As a matter of fact, we show in Section ?? of the supplementary that IPOT is a MM algorithm.
Interestingly. this scaling matrix
The idea is to re-write the OT objective as ((CTT) + (T, ab
AD (T, ab ") and upper bound the concave term by its tangent. This further supports the interest of MM for OT and UOT, and highlight an important feature of one of our contributions: designing the first Sinkhorn-like multiplicative algorithm for UOT that can be applied when the OT plan is not entropy-regularizzd.
More efficient solvers. Despite the positive experimental results of Section 4 multiplicative and regularization path algorithms can be slow, especially for large values of A. Various accelerations can be envisaged. Regarding path algorithms, the approach Maiali and Yu (2012) can compute a regularization path with precision in o(1/c) iterations. This would lead in our setting to a full complexity of (mm/e) that is even interesting approximate balanced OT. Another way to speed up computations is to use screeming. In sparse regression, this consists of eliminating during optimization components that will not belong the support of the solutions thanks to safe screening tests. Methods such as (El Ghaoui et al. 20 Wang et al. 2015 Dantas et al. 202I] can readily be adapted to our l2 or KL-penalized UOT algorithms. Finally, an other line of improvement is to consider stochastic optimization methods such as (Defazio et al., (014)). Given the particular structure of H. the complexity of stochastic updates shall be small and can lead to very efficient implementations (Nesterov] 2014).
General case and entropy-rreullaiized UOT. Following (Frogner et al. 2015 Chizat et al. Séjourné et al. 2019). general regularized UOT can be expressed as:
2018
RUOT b)
min T20
(C,T)
(TT a) + A2 A2 I,, + A D (T. ab
(13)
As it turns out, this general problem involving different regularization weights (A1, A2, Areg) can easily be addressed in our framework as well using two simple tricks. The first one consists of absorbing the regularization weights into the divergences. Indeed, many divergences are homogeneous, i.e., satisfy a relation of the form (xly) = D (>"xx Xyy) where a is divergence-specifie This holds in particular for the KL divergence (a 1) and the squared lg norm (a 1/2). The second one consists of complementing H and y with suitable terms to account for the regularization term. In the end, we may re-write Eq. (3)) into Eq. (3) with À = 1, H = [>g H >g H Arg IJT and [Ata? Ag Agg vec )T. In particular, obtain the following multiplicative update in the case of entropy-regee regularieg KL- penalize UOT:
((++1) = -diag
a
(T(A))
Art Am
diag
b T
2a T I,
(14)
T (b)
1
where K = (ab
+ X A This multiplicative update is slightly more complex than the Sinkhorn algorithms of Frogner et all (2015); Chizat et al. (2018) and as such, it might have limited practical interest but is conceptually interesting and novel. Note that balanced UOT as of Eq. (2) is simply obtained with Aveg = 0.
exp
and Aan = > A2
Aun
Non-linear UOT. Finally, we discuss how our proposed reformulation of UOT can accommodate non-linear variants in which the linear term (C,T) is replaced by htapiy//robbusteesssrooiiii

---

term, leading to problems of the form
NLUOT b) = min T20
+ a) (T In.b)
(15)
where g(-) is usually concave function, see, e.g., (Candes et al. 2008 Gasso et al. 2009). Our MM
setting can readily accommodate such a formulation by majorizing the concave terms by their tangent. The non-linearity may improve robustness w.r.t outliers and better model realistic OT problems. For instance, in real life, the costs of transporting some goods between two places can nonlinear due to economies of scale.
Broad and potential negative societal impact. The contributions in this paper are methodological and focus on a reformulation of a fundamental OT problem and adapting existing algorithms to solve In this sense, we bring more efficient solvers that run on GPU but this computational advantage can be counterbalanced by the possibillity that it brings to be applied on larger datasets. 'The application of OT in domain adaptation has shown that can be used to infer labels on samples/individuals when no labels are available, suggesting a capacity for violating user privacy. A potential application of UOT is the case where two datasets of users acquired by different methods contain some shared users. UOT can be used here to find correspondences between the users in the two datasets and also identify unique users in each dataset (those that do not receive mass).
Conclusion
In this paper, we reformulate the UOT problem as a non-negative penalized linear regression, allowing us to propose two new classes of algorithms. We first derive multiplicative algorithms for both KL and fz-penalized UOT, providing numerical solutions that are fast and easy to implement. For the specific case of fy-penalized UOT, we provide the first regularization path algorithm that computes the whole set of solutions for all the regularization parameter values. We finally build on the extensive literature inverse problem and NMF to draw some fruitful perspectives on even more efficient algorithmic solutions or the definition of new OT problems.
Acknowledgments and Disclosure of Funding
The authors want to thank Hicham Janati for interesting discussions and providing us with the experiments of convergence for the MM algorithm in the supplemental. This work is partially funded by the French National Research Agency (ANR; grants OATMIL ANR-I7-CE2S 17-CE23-0011 RAIMO ANR-20-CHLA CHIC--021----- 1-01- MULTISCALE 0006-02, 31A Côte d' Azur ANR-19-P3IA- 31A ANITI ANR-19-P134
0022-01, E4C ANR- 18-EUR- and the European Research Council (ERC: grant FACTORY-CoG-66 1866818339 1839). Furthermore, this research was produced within the framework of Energy4Climate Interdisciplinary Center (E4C) of IP Paris and Ecole des Ponts ParisTech. This action benefited from the support of the Chair "Challenging Technology for Responsible Energy" led by I'X- Ecole Polytechnique and the Fondation de I' Ecooee Polytechnique, sponsored by TOTAL.
References
Arjovsky, M., Chintala, and L. Bottou (2017). Wasserstein Generative Adversarial Networks. In International Conference on Machine Learning, Volume 70, pp. 214-223.
Benamou, J.-D. (2003). Numerical resolution of "unbalanced" mass transport problem. ESAIM: Mathematical Modelling and Numerical Analysis 37(5), 851-868.
Blondel, M., V. Seguy, and A. Rolet (2018). Smooth and Sparse Optimal Transport. In International Conference on Artificial Intelligence and Statistics, pp. 880-889.
Bonneel, N. and D. Coeurjolly (2019). SPOT: Sliced Partial Optimal Transport. ACM Transactions on Graphics (SIGGRAPH) 38(4).
Byrd, R. H., Lu, Nocedal, and C. Zhu (1995). Limited Memory Algorithm for Bound Constrained Optimization. SIAM Journal Scientific Computing 16(5), 1190-1208.
10

---

Caffarelli, L. A. and R. J. McCann (2010). Free boundaries in Optimal Transport and Monge-Ampère obstacle problems. Annals of Mathematics (2)( 673-730.
Candes, E. J. M. B. Wakin, and S. P. Boyd (2008). Enhancing Sparsity by Reweighted & Minimiza- tion. Journal of Fourier analysis and applications 14(5-6), 877-905.
Chizat, L., G. Peyré, B. Schmitzer, and F.-X. Vialard (2018). Scaling algorithms for Unbalanced Optimal Transport problems. Mathematics of Computation 87(314), 2563-2609.
Courty, N., R. Flamary, D. Tuia, and A. Rakotomamoniy (2017). Optimal Transport for Domain Adaptation. IEEE Transactions on Puttern Analysis and Machine Intelligence 39(9), 1853-1865.
Cuturi, M. (2013). Sinkhorn distances: Lightspeed computation of Optimal Transport. Advances in Neural Information Processing Sytems 26, 2292-2300.
Dantas, C. E. Soubies, and C. Févotte (2021). Safe Screening for Sparse Regression with the Kullback-Leibler Divergence. In IEEE International Conference Acoustics, Speech and Signal Processing, pp. 5544-5548.
De Pierro, A. R. (1993). On the relation between the ISRA and the EM algorithm for positron emission tomography. IEEE transactions on Medical Imaging 12(2). 328-333.
Defazio, A., F. Bach, and Lacoste-Julien (2014). SAGA: fast incremental gradient method with support for non-strongly convex composite objectives. In Advances in Neural Information Processing Systems.
Dhillon, I. S. and S. Sra (2005). Generalized nonnegative matrix approximations with Bregman divergences. In Advances in Neural Information Processing Systems, Volume 18.
Efron, B., T. Hastie, I. Johnstone, and Tibshirani (2004). Least Angle Regression. Annals of statistics 32(2), 407-499.
El Ghaoui, L., V. Viallon, and T. Rabbani (2012). Safe Feature Elimination for the LASSO and Sparse Supervised Learning Problems. Pacific Journal of Optimization 866776988.
Févotte, C. and J. Idier (2011). Algorithms for nonnegative matrix factorization with the 3-divergence. Neural computation 23(9), 2421-2456.
Figalli, A. (2010). The Optimal Partial Transport Problem. Archive for Rational Mechanics and Analysis 195(2). 533-560.
Flamary, R., N. Courty, A. Gramfort, M. Z. Alaya, A. Boisbunon, S. Chambon, L. Chapel, A. Coren- flos, K. Fatras, N. Fournier, et al. (2021). POT: Python optimal transport. Journal of Machine Learning Research 22(78), 1-8.
Frogner, C., C. Zhang, H. Mobahi, M. Araya, and T. A. Poggio (2015). Learning with a Wasserstein Loss. In Advances in Neural Information Processing System, pp. 2053-2061.
Gasso, G., A. Rakotomamoniy, and S. Canu (2009). Recovering sparse signals with a certain family of nonconvex penalties and DC programming. IEEE Transactions Signal Processing 57(12), 4686-4698.
Hastie, T., S. Rosset, R. Tibshirani, and J. Zhu (2004). The entire regularization path for the Support Vector Machine. Journal of Machine Learning Research 1391-1415.
Ho, N., X. L. Nguyen, M. Yurochkin, H. H. Bui, V. Huynh, and D. Phung (2017). Multilevel Clustering via Wasserstein Means. In International Conference Machine Learning, Volume 70, pp. 1501-1509.
Hunter, D. R. and K. Lange (2004). tutorial on MM algorithms. The American Statistician 58(1), 30-37.
Janati, H. (2021). Advances in Optimal transport and applications to neuroscience. Ph. D. thesis,
Institut Polytechnique de Paris.
11

---

Janati, H.. B. Muzellec, G. Peyré, and M. Cuturi (2020). Entropic optimal transport between unbalanced gaussian measures has a closed form. Advances in Neural Information Processing Systems 33.
Kantorovich, L. (1942), On the transfer of masses (in Russian). Doklady Akademii Nauk 2, 227-229.
Kusner, M., Y. Sun, N. Kolkin, and K. Weinberger (2015). From word embeddings to document distances. In International Conference on Machine Learning, pp. 957-966.
LeCun, Y., C. Cortes, and C. Burges (2010). MNIST handwritten digit database. ATT Labs [Online]. Available: http://yann. www.n. wwccn.cco//xxd/mmiss..
Lee, D. and Seung (1999). Learning the parts of objects by non-negative matrix factorization. Nature (67(5),,), 788-791.
Lee, D. and H. Seung (2001). Algorithms for Non-negative Matrix Factorization. In Advances in Neural Informaiion Processing Sytems, Volume 13.
Liero, M., A. Mielke, and G. Savaré (2018). Optimal entropy- problems and a new Hellinger-
Kantorovich distance between positive measures. Inventiones mathematicae 211(3), 969-1117.
Mairal, and B. Yu (2012). Complexity Analysis of the Lasso Regularization Path. International Conference on Machine Learning, pp. 1835-1842.
Maretic, H. P.. M. E. Gheche, G. Chierchia, and P. Frossard (2019). GOT: An Optimal Transport framework for Graph comparison. In Advances In Neural Information Processing Systems, Volume 32.
Massias, M., A. Gramfort, and J. Salmon (2018). Celer: Fast Solver for the Lasso with Dual Extrapolation. In International Conference on Machine Learning, Volume 80, pp. 3321-3330.
Mukherjee, D., A. Guha, J. Solomon, Y. Sun, and M. Yurochkin (2020). Outlier-Robuu Optimal Transport. Technical report, arXiv preprint arXiv:2012.07363.
Nesterov, Y. (2014). Subgradient methods for huge-scale optimization problems. Mathematical Programming 146(1), 275-297.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Pretten- hofer, R. Weiss, V. Dubourg, et al. (2011). Scikit-learn: Machine learning in Python. Journal of Machine Learning Research 12, 2825-2830.
Rabin, J., S. Ferradans, and N. Papadakis (2014). Adaptive color transfer with relaxed optimal transport. In IEEE International Conference on Image Processing, pp. 4852-4856.
Richardson, W. H. (1972). Bayesian-based iterative method of image restoration. JoSA 62(1), 55-59.
Sato, R., M. Yamada, and H. Kashima (2020). Fast Unbalanced Optimal Transport Tree. In Advances in Neural Information Processing Systems, Volume 33.
Séjourné, T., J. Feydy, F.-X. Vialard, A. Trouvé, and G. Peyré (2019). Sinkhorn divergences for Unbalanced Optimal Transport. arXiv preprint arXiv:1910.122 12958.
Sinkhorn, R. and P. Knopp (1967). Concerning nonnegative matrices and doubly stochastic matrices. Pacific Journal of Mathematics 2/(2), 343-348
Sun, Y., P. Babu, and D. Palomar (2017). Majorization-mmniiizzttoo algorithms in signal processing, communications, and machine learning. IEEE Transactions on Signal Processing 65(3), 794-8194 816..
Thibault, A., L. Chizat, C. Dossal, and N. Papadakis (2021). Overrelaxed Sinkhorn-Knopp Algorithm for Regularized Optimal Transport. Algorithms 14(5), 143.
Vayer, T., Chapel, R. Flamary, R. Tavenard, and N. Courty (2019). Optimal transport for structured data with application on graphs. In International Conference on Machine Learning, pp.
Villani, C. (2009). Optimal Transport: Old and New, Volume 338. Springer Berlin Heidelberg.
12

---

Virtanen, P. R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau. E. Burovski,
P. Peterson, W. Weekesser, J. Bright, S. van der Walt. M. Brett, Wilson, K. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C. J. Carey, I. Polat, Y. Feng, E. W. Moore, J. VanderPlas, D. Laxalde, Perktold, Cimrman, I. Henriksen, E. A. Quintero, C. R. Harris, A. M. Archibald, A. H. Ribeiro, F. Pedregosa, P. van Mulbregt, and SciPy 1.0 Contributors (2020). SciPy 1.0: Fundamental Algorithms for Scientific Computing Python. Nature Methods 17, 261-272.
Wang. J., P. Wonka, and J. Ye (2015). Lasso screening rules via dual polytope projection. Journal of Machine Learning Research 16(1), 1063-1101.
Xiao, H., K. Rasul, and R. Vollgraf (2017). Fashion-MNIST: a Novel Image Dataset for Benchmarking Machine Learning Algorithms.
Xie, Y., X. Wang, Wang. and H. Zha (2020). A fast proximal point method for computing exact Wasserstein distance. In Uncertainty in Artificial Intelligence, pp. 433-453.
Yang, Z. and E. Oja (2011). Unified Development of Multiplicatiive Algorithms for Linear and Quadratic Nonnegative Matrix Factorization. IEEE Transactions on Neural Networks 22, 1878 1891.
13


Reviews:
[{"review_id": "ewawbkMea0J", "reviewer": "Reviewer_7tpv", "summary": "In this work, the authors propose to phrase the Unbalanced Optimal Transport problem under general Bregman divergence penalty as a non-negative penalized linear regression problem. They exploit an existing algorithm in Non-negative Matrix Factorization problem (NMF), namely Majorization-Minimization, in order to solve the obtained problem. In the case of KL-penalized UOT, the obtained algorithm is a variation of the Sinkhorn algorithm. However, this algorithm solves a different algorithm than Sinkhorn, namely KL-penalized UOT as opposed to OT with KL penalized transport plan.  For the l2-penalized UOT, the authors design an algorithm, based on LARS, that solves the problem for all values of lambda, from 0 to infinity.\n\nThe running time of different variants of their proposed algorithm are compared on a toy problem versus the current state-of-the-art method that uses BFGS. Overall, for l2-penalized UOT, the proposed method is significantly faster than BFGS. For KL-penalized UOT, it depends on the value of lambda.\n", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "pros:\n- The connection made between UOT and non-negative linear regression allows for the use of well-studied efficient algorithms.\nFor KL-penalized UOT, the algorithm solves the problem without needing to regularize the transport plan, as opposed to using Sinkhorn divergence for UOT as done in [1].\n- For l2-regularized UOT, the proposed algorithm can be solved once for all values of the regularization parameter, and is shown to be faster than the state-of-the-art method that uses BFGS.\n\ncons:\n- The improved performance of the algorithm is shown only in the case of l2-penalized UOT, whose practical interest is not clear compared to KL-penalized UOT, which is more commonly seen in the literature.\n- In the case of KL-penalized UOT, I expected a deeper analysis of the connection with the Sinkhorn algorithm and the Sinkhorn divergence based UOT of [1], e.g., how similar do they become for asymptotic values of lambda and epsilon (epsilon being the reg. parameter for the entropy of the transport plan).\n\nTypos:\n\nl.319 balanced OT \n\n[1] Sinkhorn Divergences for Unbalanced Optimal Transport, Séjourné et al.\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors adequately addressed the limitations and potential negative societal impact.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3h"}, "extra_text": {}}, {"review_id": "SkL-orYXG9x", "reviewer": "Reviewer_42ef", "summary": "This paper proposes a new method for Unbalanced Optimal Transport (UOT) without a regularization term. The solution consists of two steps; (i) formulate the UOT as a non-negative penalized linear regression problem (ii) solve it using the Majorization-Minimization (MM) algorithm. The authors derive specific algorithms for several penalty functions. Furthermore, a regularization path tracking method based on the relationship between the proposed formulation and lasso is proposed. The authors validate the effectiveness of the proposed method by numerical experiments using synthetic datasets and the MNIST dataset.", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "### Strengths\n\n* The idea of formulating UOT as a non-negative penalized linear regression problem and applying the MM algorithm is very interesting.\n\n* Also, there is originality in the regularized path tracking algorithm based on the above formulation.\n\n* This paper is well-written; the problems, the backgrounds, and the contributions are clearly stated with appropriate citations. \n\n### Weaknesses\n\n* The analysis of the time complexity is insufficient both theoretically and experimentally. The theoretical analysis of the time complexity for each iteration and convergence is not sufficiently done for any of the methods, and the superiority of the proposed method cannot be judged. As for the experiments, the scalability of the proposed method cannot be judged from the experimental results because the comparison with the existing methods is made only when \\lambda is changed, and there is no comparison when n or m is changed. Although the authors do not claim that the computation time is an advantage of the proposed method, these analyses are necessary to understand the characteristics of the proposed method. \n\n* In the introduction, the authors state that one of the advantages of the proposed method is fast computation by GPUs, but this claim has not been verified theoretically or experimentally.\n\n* In the introduction, the ease of implementation is mentioned as an advantage of the proposed method, but since the L-BFGS method can also be easily implemented by using scipy, there seems to be no significant difference.\n\n### Other comments\n\n* Regularized path in OT is interesting, but it is a bit unclear how it can be useful. I would like to know if there are any specific use cases.\n\n* In the experiments, the authors compare their method with the BFGS method in Blondel et al. (2018).  Blondel et al. (2018) propose two methods in the paper, one for solving the dual and one for solving the semi-dual.  Which method is used in the experiment? In addition, Blondel et al. (2018) also describe an optimization method using alternating minimization. Because this method resembles the proposed method, shouldn't it be compared with this method as well?\n\n* Section 5 contains various speed-up proposals for the proposed algorithm.  A thorough theoretical and experimental evaluation of these proposals will make this paper more interesting.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The limitation of this paper is that it has little advantage over existing methods in terms of computation time, as the authors also state in the checklist.  The societal impacts are fully discussed at the end of Section 5.", "review_text": "", "extra_scores": {"ethical_concerns": "N/A", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "7"}, "extra_text": {}}, {"review_id": "IhqkJIDyyYF", "reviewer": "Reviewer_1XGy", "summary": "Authors in this paper flattened the matrix-form unbalanced OT and reformulated it as a linear regression regularized by either the KL divergence or $\\ell_2$ between the transport plan and two marginals. Then the problem converts to an MM problem and the authors derived the optimization steps for solving them, especially the one with $\\ell_2$ in which case the authors proposed an algorithm to regress the optimal combinations of the transport plan and the regularization strength. Numerical results show that the new formulation leads to faster computation.", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "+\nThe proposed method is novel and all the claims in the paper seem well supported by theoretical and empirical results.\n\n-\nReformulating $\\ell_2$-regularized UOT to a Lasso regression inevitably promotes sparsity. Adding the transportation cost $c^Tt$, the resulting transporting map $t^k$ concentrates on the nearest samples of the two domains like what $Figure\\ 2$ shows, especially when $\\lambda$ is small. When $\\lambda_1 = \\lambda_2$, it becomes more apparent that $Algorithm\\ 1$ is selecting a subset of the samples and finding a balanced transport. Although it still aligns with the original formulation of unbalanced OT (Benamou, 2003), but my understanding of the practicality of unbalanced OT is that it still transports all the mass but allows the mass to grow or shrink during transport. Now, the proposed method seems more like to transport a subset of mass controlled by $\\lambda$. \n\nI am also reluctant to accept the wording that $\\lambda \\rightarrow \\infty$ recovers a balanced OT because balanced OT requires the mass on both sides to equal, i.e. $ \\mathbb{1}a =  \\mathbb{1}b$ and such a T does not exist when $\\lambda \\rightarrow \\infty$ and $ \\mathbb{1}a \\neq  \\mathbb{1}b$. In this sense, the experiment in Figure 2 feels cheesy to me because the two marginals are balanced, and $Algorithm\\ 1$ feels more like gradually finding a balanced OT by including more samples and excluding fewer outliers from the two domains. \n\nThe motivation behind computing the \"regularization path\" is not clear in the beginning of the paper. Only when I saw the application to unbalanced DA did I get to know its usage but I'm still not sure what motivated the authors to get such as path. It seems the paper has a very similar motivation to those of Mukherjee et al. (line 416) and Balaji et al.\n\nBalaji, Yogesh and Chellappa, Rama and Feizi, Soheil, Robust Optimal Transport with Applications in Generative Modeling and Domain Adaptation, Advances in Neural Information Processing Systems, 2020\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3"}, "extra_text": {}}, {"review_id": "IAteZO5Gzy6", "reviewer": "Reviewer_9N5u", "summary": "The authors propose to reformulate the unbalanced Optimal Transport (e.g., Kullback-Leibler and $\\ell_2$-regularization) and non-negative penalized linear regression with $\\ell_1$-regularization. For $\\ell_2$-regularized UOT, the authors derive the regularization path.\n", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "It is interesting to reformulate the UOT as a non-negative penalized linear regression with $\\ell_1$ regularization and leverage the tools in Lasso problem to derive the regularization path.\n\nThe majorization-minimization (MM) algorithmic approach for UOT seems elegant with the multiplicative update (naturally handing the non-negative constraint). However, there is no theoretical analysis to compare this MM algorithmic approach for the proposed reformulation with approaches in the literature (for KL/$\\ell_2^2$ regularized UOT). Although the authors did the empirical comparison in Figure 3 (and recall that the complexity of L-BFGS is about $O(n^2)$), it seems unclear why the proposed MM algorithmic approach gives more advantages. It is better in case the authors elaborate with more details.\n\nThe regularization path for $\\ell_2$-penalized UOT is interesting. However, it seems that the authors may need to discuss Algorithm 1 in more detail. For examples,\n+ Can the stoping condition in the while-loop be obtained?\n+ How to update the set A? (How to compute the set $A_{k+1}$?\n\nFor the right plot in Figure 1 (left), it may mislead the readers for the curve between $\\lambda_7$ and $\\infty$? (monotonically increasing?)\n+ Do the authors compute the standard OT at $\\infty$ or use Algorithm 1?\n\nIn line 290-292, the authors discuss the contribution with the SInkhorn-like multiplicative algorithm for UOT? However, it is unclear why the Sinkhorn-like algorithm is important and favoring? What is its advantage to compare with the L-BFGS (e.g., in Blondel et al.)\n\nOverall, I think the submission is interesting. However, there are some points (as mentioned above) that may be needed further elaboration. Besides the regularization path results (for $\\ell_2$ regularization), it is still unclear yet about the advantage of the MM algorithmic approach for the proposed reformulation with approaches in the literature.\n\n---------After the rebuttal-----------\n\nI thank the authors for the response. Overall, I think the proposed reformulated UOT to connect with the LASSO problem is interesting. I hope that the authors address those raised concerns, especially those about the MM algorithm in the updated version. So, I keep learning on the positive side. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes", "review_text": "", "extra_scores": {"ethical_concerns": "N/A", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "8 hours"}, "extra_text": {}}]

Author rebuttal:

Meta-review: The authors consider the unbalanced optimal transport (OT) problem and its semi-relaxed variants. They then consider various optimization approaches that come from the mirror descent family by choosing an appropriate Bregman divergence. In particular, they consider the entropic and the ell2 mirror maps. 

In the ell2 case, the authors connect the optimization problem to the LARS optimization, which has been quite successful in the sparse recovery literature. They propose a regularization path or homotopy algorithm for the unbalanced OT plans and demonstrate its performance numerically.

The developments in the paper mostly rely on recognizing the LASSO structure in the unbalanced OT problem and then applying existing techniques to it. In the worst case, this approach has exponential complexity but the authors argue that the [Mairal and Yu 2012] approach can be applied to obtain eps-approximate primal dual gap as an afterthought. Unfortunately, there is also no follow up to back up how this epsilon guarantee to the quality of the OT map. 




Decision: Accept (Poster)