All You Need is Resistance: On the Equivalence of Effective Resistance and Certain Optimal Transport Problems on Graphs
Sawyer Robertson Zhengchao Wan Alexander Cloninger


	
	Estimation Network Design framework for efficient distributed optimization

    Mattia Bianchi and  Sergio Grammatico
     
	M. Bianchi is with the Automatic Control Laboratory (IfA), ETH Zürich, Switzerland (). S. Grammatico is with the Delft Center for Systems and Control (DCSC), TU Delft, The Netherlands (). This work is supported by NWO under research project OMEGA (613.001.702), by the ERC under research project COSMOS (802348) and by ETH Zürich funds.

    Received 28 February 2024 / Accepted 23 April 2024
==================================================================================================================================================================================================================================================================================================================================================================================================



    
        The fields of effective resistance and optimal transport on graphs are filled with rich connections to combinatorics, geometry, machine learning, and beyond. In this article we put forth a bold claim: that the two fields should be understood as one and the same, up to a choice of p. We make this claim precise by introducing the parameterized family of p-Beckmann distances for probability measures on graphs and relate them sharply to certain Wasserstein distances. Then, we break open a suite of results including explicit connections to optimal stopping times and random walks on graphs, graph Sobolev spaces, and a Benamou-Brenier type formula for 2-Beckmann distance. We further explore empirical implications in the world of unsupervised learning for graph data and propose further study of the usage of these metrics where Wasserstein distance may produce computational bottlenecks.
    

    
    Keywords: effective resistance, optimal transport, graph Laplacians, convex optimization, random walks
    

    
    MSC: 65K10, 05C21, 90C25, 68R10, 05C50
    
	

	

§ INTRODUCTION AND NOTATION

	
	

 §.§ Introduction

	
	The theory of optimal transportation (OT), which traces its roots to the Monge formulation as early as 1781 <cit.> as well as the Kantorovich relaxation articulated in 1942 <cit.>, is typically framed as an optimization problem for some convex or affine function which measures the cost of mass transportation over all feasible plans which provide recipes for transporting some mass distribution to another. OT has received significant attention due to its connections to geometry <cit.>, partial differential equations <cit.>, applied mathematics <cit.>, and many other fields. For discrete measures and domains in particular, OT and Wasserstein distance have been used and studied for a variety of purposes including document retrieval <cit.>, statistics <cit.>, image registration <cit.>, distance approximation <cit.>, political redistricting <cit.>, graph neural networks <cit.>, and graph Ollivier-Ricci curvature <cit.>.
	
	When the underlying metric space of the transportation process is a graph, much attention has been focused on the 1-Wasserstein metric, which can be informally set up as a linear program of the following form:
		
    𝒲_1(α,β) = inf{∑_ijπ_ij d(i, j) : π∈Π(α,β)},

	where Π(α,β) consists of the couplings of two generic probability distributions α,β on the nodes of a graph and d(i, j) is the shortest path metric on the vertex set V (see <ref> for a complete description). It so happens on graphs that this metric can also be written in a fashion often termed the Beckmann formulation:
		
    𝒲_1(α,β) = inf{∑_e |J(e)|w_e : J:E→ℝ, B J = α-β}

	where B is the graph incidence matrix, J is any map on the edges of the graph and which is identified as a vector in ℝ^|E|, and w_e is the nonnegative weight for edge e∈ E (see <ref> for a complete description). This formulation is named after Martin Beckmann who studied dynamic formulations of OT beginning in the 1950s <cit.>, and whose models for OT in the continuous setting <ref> closely mirrors in a discrete sense. The variable J is sometimes called a feasbile flow between α and β (see <ref>), and as such <ref> is also sometimes deemed a flow-based formulation for OT on graphs.
	
	Seemingly separately, effective resistance (ER) on graphs has long been studied for its metric properties between nodes and its connections to many other problems in graph theory. Recall that for two nodes i, j∈ V, and corresponding standard basis vectors δ_i, δ_j∈ℝ^|V|, the ER between them is given by r_ij = (δ_i - δ_j)^T L^† (δ_i - δ_j) where L^† is the Moore-Penrose psuedo inverse of the graph Laplacian. Applications and research on ER include foundational work by Spielman and Srivastava on spectral graph sparsifcation using ER <cit.>, as well as spectral clustering models leveraging computational advantages of ER solvers <cit.>. ER has also been used for computing and studying statistics associated to random walks on graphs <cit.>. Interestingly, in the case of directed and connection graphs, the opposite is the case as random walk statistics have been used to obtain versions of ER <cit.>. Moreover, there exist geometric applications involving ER as an ingredient for defining discrete Ricci curvature <cit.>. For other references which study and utilize effective resistance in various contexts, see, e.g., <cit.>.
	
	This paper is built on an observation which seems to have not yet appeared in the field: that effective resistance, when extended to probability measures on a graph G, can be realized as an OT distance obtained from a Beckmann-type problem with a 2-norm penalty. The resulting metric, while removed from a conventional coupling-based formulation of OT, is intricately connected to the theory of random walks on graphs, the graph Laplacian matrix, and graph Sobolev spaces. 
		Let G=(V,E, w) be a weighted connected graph with Laplacian matrix L and psuedoinverse L^†. Denote the oriented vertex-edge incidence matrix of G by B, and let α,β be fixed probability measures on V regarded as vectors in ℝ^|V|. For 1≤ p<∞, define ℬ_p(α,β) by the following convex optimization problem:
            
    ℬ_p(α,β)^p = inf{∑_e∈ E|J(e)|^p w_e : J:E→ℝ, BJ = α-β}.

        Then it holds:
            
                
  * When p=1, ℬ_1(α,β) = 𝒲_1(α,β).
                
  * When p=2, ℬ_2(α,β) = (α-β)^T L^† (α -β).
            
	
	
	The metric introduced in <ref> for p=2 is a (squared) norm distance, so its properties and structures are fundamentally different from the (nonlinear) Wasserstein metric in general. However, this metric seems to posses a rich array of properties and characterizations which demonstrate strong homophily with the classical results in OT. This metric also shows promise for use in kernel-based learning methods for data defined on graphs, and in terms of complexity, is significantly less expensive to compute at scale compared to the Wasserstein metric.
	
	Our main contributions are summarized as follows. 
	
	
		
  * We introduce the notion of p-Beckmann distance between probability measures, and derive various bounds between p-Beckmann distances and p-Wasserstein distances for measures on general graphs (e.g., <ref>).
        
  * We provide an explicit formula for ℬ_p on trees (<ref>).
		
  * We obtain a generalized commute time formula for 2-Beckmann distance, or measure effective resistance, in terms of optimal access times of the simple random walk on a graph (<ref>).
		
  * We realize 2-Beckmann distance as a negative Sobolev-type semi-norm distance and use this property to obtain a Benamou-Brenier type formula for 2-Beckmann distance (<ref>).
		
  * We apply 2-Beckmann distance to an unsupervised kernel-based clustering method involving handwritten digit data, and compare the results with the 2-Wasserstein kernel; and thereby establish empirical evidence for the usage of a 2-Beckmann kernel on a drop-in basis in place of a 2-Wasserstein kernel for datasets defined on graphs (<ref>).
	

	

 §.§ Outline of this paper


	In <ref>, we provide relevant mathematical background and notation.
	
	In <ref>, we delve into the general properties of the p-Beckmann problem, including duality theory in <ref> and estimates between Beckmann and Wasserstein distances <ref>. In <ref>, we highlight some specialized results in the cases of paths and trees.
	
	In <ref> we focus exclusively in the properties and theory of the 2-Beckmann distance. In <ref> we frame the metric as effective resistance between measures and relate it to optimal stopping times and access times between probability measures on G. In <ref> we frame the metric as a negative Sobolev-type seminorm and obtain a graph Bemaou-Brenier-type formula. Finally, in <ref> we apply the 2-Beckmann distance to an unsupervised learning problem and explore various empirical observations and implications. 
	
	Lastly, in <ref> we provide proofs of the results mentioned in the main component of the paper.
	
    

 §.§ Related work

 
 	In this section we wish to highlight some particularly relevant works on which this paper is based in parts, and in doing so, better situate the novelty of our contributions in a broader context.
 	
 	The most relevant investigation of OT on graphs, including regularization and convex duality, appears in work by Essid and Solomon <cit.> and further studies in the field can be found in <cit.>.
 	
 	The introduction of p-Beckmann distances in <ref> and its study as an effective resistance metric in <ref> should be understood as extensions of the notion of p-resistances from work of Alamgir and Luxburg <cit.> to general probability measures (from the case of nodes, i.e., Dirac measures). Such generalizations of resistances and their implications in unsupervised models were also considered by Nguyen and Hamitsuka <cit.>, and Saito and Herbster <cit.> but only between nodes. The notions of measure access times and optimal stopping times have been studied at length by Lovász, Winkler and Beveridge across several papers <cit.>.
 	
 	Separately, the notion of using negative Sobolev distance as a linearization of quadratic Wasserstein distance was studied by Greengard, Hoskins, Marshall, and Singer in <cit.> and is a homophilous framework appearing in the continuous setting. The proof of the Benamou-Brenier-type formula for 2-Beckmann distance is a direct discretization of the proof appearing in <cit.>. It is also notable that our model of 2-Beckmann distance as a negative Sobolev distance on graphs is parallel to the works by Le, T. Nguyen, Phung, V. A. Nguyen, and Fukumizo and appearing in <cit.>; but our definitions differ slightly from theirs. Linearized OT and its usage in the unsupervised setting is also explored by Moosmüller and AC <cit.>, and results therein form a basis for our <ref> on linear separation of measures in the 2-Beckmann metric space.
 
	

 §.§ Mathematical background

	
	Let G=(V, E, w) be a graph, where V = {1, 2, …, n} is the set of vertices, E ⊂V 2 is a set of undirected edges of cardinality m≥ 0, and w = (w_ij)_i, j∈ V is a choice of real edge weights satisfying w_ij≥ 0, w_ij = w_ji, and w_ij > 0 if and only if {i, j}∈ E. In order to ensure the feasibility of optimal transportation problems and simplify the exposition in places, we assume that G is finite, has no multiple edges or loops, and is connected. 

	For our purposes, a path in G is an ordered sequence of nodes P = (i_0, i_1,…, i_k) such that i_ℓ∼ i_ℓ+1 for 0≤ℓ≤ k-1. If P is a path, |P| is the (integral) length of the path, i.e., |P|=k, and |P|_w is the weighted length of P, or the sum of the weights of edges it contains; i.e. |P|_w = ∑_ℓ=0^k-1w_i_ℓ i_ℓ+1. For any two nodes i, j∈ V, 𝒫(i, j) is the set of all paths which begin at i and end at j. 
	
	In this paper, we will consider several different metrics on V. A metric on V is a map κ:V× V→ℝ such that 
		
			
  * (Nonnegativity) κ(i, j)≥ 0 for each i, j∈ V,
            
  * (Definiteness) κ(i, j)= 0 if and only if i=j,
			
  * (Symmetry) κ(i, j) = κ(j, i) for each i, j∈ V,
			
  * (Triangle inequality) κ(i, k)≤κ(i, j) + κ(j, k) for each i, j, k∈ V.
		 
	We will primarily use variants of the shortest path metric on V and the metrics induced by effective resistance.
	
	Let 1≤ p <∞. Define the weighted p-shortest path metric on V by the formula
			
    d_p(i, j) = inf{(∑_ℓ=0 ^k-1w_i_ℓ i_ℓ+1^p)^1/p : P=(i=i_0, i_1,…, i_k=j)∈𝒫(i, j)}

	

	It is relatively straightforward to verify that d_p satisfies the properties necessary for a metric on V for each p, w. 
 
    We define the Adjacency matrix A∈ℝ^n× n entrywise by
		
    A_ij   =
    				w_ij    if i∼ j
    
    				0     otherwise.

	where the notation i∼ j means {i, j}∈ E. For each i∈ V we define its degree d_i = ∑_j∼ i w_ij. We will also use the diagonal degree matrix D = diag(d_1,…, d_n)∈ℝ^n× n and the diagonal edge weight matrix W = diag(w_e_1, …, w_e_m)∈ℝ^m× m.
	
	For technical purposes, we may ocassionally refer to the index-oriented and bi-oriented edge sets E' and E”, respectively, defined below:
		
    E'    = {(i, j) : i, j∈ V, i∼ j, i < j}
    
    			E”   = {(i, j), (j, i): i, j∈ V, i∼ j}.

	
	We define the incidence matrix B∈ℝ^n× m, with rows indexed by V and columns indexed by E', by the formula:
		
    B_i, e_j   =
    				1     if e_j = (i, ·)
    
    				-1     if e_j = (·, i)
    
    				0    otherwise.

	We also occasionally use the matrix B∈ℝ^n× 2m,  with rows indexed by V and columns indexed by E”, defined by an identical formula as in <ref>. We define the Laplacian matrix L by the formula L = D-A, or L=BWB^T. Note that BWB^T = 2L, where W∈ℝ^2m× 2m is the diagonal matrix of edge weights of E”. We recall that L is symmetric, positive semi-definite, and has a set of nonnegative eigenvalues
		
    0=λ_0≤λ_1≤…≤λ_n.

	We use Λ = diag(λ_0,λ_1,…,λ_n)∈ℝ^n× n for the diagonal matrix of eigenvalues, and U = [ u_1 u_2   ⋯ u_n ] for the respective orthonormal eigenvectors of L. The spectral decomposition of L is given by the equation L=UΛ U^T.
	
	For any matrix X∈ℝ^ℓ× k, we will use the notation X^†∈ℝ^k×ℓ for its Moore-Penrose psuedoinverse. We recall that X^† is the unique matrix satisfying the following four properties <cit.>:
		
			
  * XX^† X = X
			
  * X^† XX^† = X^†
			
  * (XX^†)^T = XX^†
			
  * (X^† X)^T = X^† X
		
	where (·)^T indicates matrix transpose.
	
	In the case of the Laplacian matrix L, we have 
		
    L^† = U Λ^-1 U^T

	where, with an abuse of notation, we write
		
    Λ^-1 = diag(λ_0 = 0,λ_1^-1, …,λ_n^-1)∈ℝ^n× n.

	Similarly, we write L^-1/2 = U Λ^-1/2 U^T.
	
	For a finite set S, we define ℓ_2(S) to be the linear space of functions f:S→ℝ with the standard Euclidean inner product ⟨ f, g⟩_ℓ_2(S) = g^T f. We identify ℓ_2(V) with ℝ^n and ℓ_2(E') with ℝ^m, respectively.
	
	For a vector x∈ℝ^n and 1 ≤ p <∞, we define the vector p-norm by
		
    x_p = (∑_i=1^n |x_i|^p)^1/p.

	For any finite set S, in a similar manner as before, we can define ℓ_p(S) to be the linear linear space of functions f:S→ℝ with the norm ·_p. If in some cases we wish to refer simply to the linear space of functions on S, such as when we employ multiple norms at once, we use the notation ℓ(S) and specify the norm or inner product in context.
	
	If J∈ℓ(E'), then we can weight the norm ·_p by the edge weights as follows:
		
    J_w, p   =(∑_e∈ E'|J(e)|^p w_e)^1/p,

	and similarly for ℓ(E”). We will use the subscript of w to indicate the presence of edge weights.
	
	As a matter of notation, functions in ℓ(V) will be denoted with lowercase Latin and Greek letters (in the case of measures), and functions in ℓ(E') and ℓ(E”) will be denoted with uppercase Latin letters.

	Note that the matrices B, B^T, and L act on their respective function spaces by matrix multiplication. We summarize their pointwise formulas below.
		
    (BJ)(i)    = ∑_e∈ E': e=(i,·) J(e) - ∑_e∈ E': e=(·, i) J(e),     i∈ V, J∈ℓ(E')
    
    			(B^Tf)(e=(i, j))    = f(i) - f(j),     e∈ E', f∈ℓ(V)
    
    			(Lf)(i)    =∑_j∼ i (f(i) - f(j))w_ij,      i∈ V, f∈ℓ(V)

	It is worth noting that in some contexts (e.g., <cit.>), B(·) and B^T(·) are considered the graph equivalents of the divergence and gradient operators, respectively. We will also occasionally use the Dirichlet energy functional for f∈ℓ(V), defined by the quadratic form of the Laplacian matrix:
		
    f^TLf =∑_(i, j)∈ E'w_ij(f(i) - f(j))^2 = B^T f_w, 2^2

	
	We define the probability measure simplex 𝒫(V) by the set
		
    𝒫(V) := {α∈ℓ(V) : α≥ 0, ∑_i∈ Vα(i) = 1}.

	For i∈ V, δ_i is the Dirac or unit measure at node i, identified with the i-th standard basis vector in ℝ^n.
	
	The primary objects of study for this paper are two different optimal transportation distances between probability measures on V. The first, which is classical and extensively studied on graphs, is p-Wasserstein distance.
	
	Let α,β∈𝒫(V), 1 ≤ p <∞, and k any metric on V. Define the set of transportation couplings between α and β, denoted Π(α,β), by the following set
			
    Π(α,β) = {π∈ℝ^n× n : π≥ 0, π1 = α,1^Tπ = β^T},

		where 1∈ℝ^n is the vector containing all ones. We define the (k,p)-Wasserstein distance between two probability measures , denoted 𝒲_k,p(α, β) by the following optimization problem:
			
    𝒲_k,p(α, β) = inf{(∑_i, j∈ Vπ_ijk(i, j)^p)^1/p :  π∈Π(α,β)}.

	
	
	As a matter of convention, when k = d_1, that is, the weighted 1-shortest-path metric, we write 𝒲_d_1,p = 𝒲_p. The second optimal mass transportation distance is much less studied in its own right except in the cases of p=1 and to a seemingly much lesser extent, p=2.

	Let 1 ≤ p <∞ and α,β∈𝒫(V). Define the set of feasible edge flows between α and β, denoted 𝒥(α,β), by the affine region
			
    𝒥(α,β)    = {J∈ℓ_2(E') : BJ=α-β}.

		Then the p-Beckmann distance between α,β, denoted ℬ_p(α,β) is given by the following constrainted norm optimization problem:
			
    ℬ_p(α,β)    = inf{J_w,p : J∈𝒥(α,β)}

	

    

	Note that 𝒥(α,β)≠∅ for any α,β since the column space of B is exactly the mean zero functions on V. The p-Wasserstein and p-Beckmann problems are essentially two different perspectives on optimal transportation for measures on graphs. The “Wasserstein philosophy" being that mass transportation is tracked between all pairs of nodes, and is accounted according to d_ij^p. We can also think of this philosophy as one of mass teleportation- the primal variable, or coupling, does not reveal how mass π_ij moves from i to j, only that it did so along a shortest path somewhere, or a combination of shortest paths. The “Beckmann philosophy" is the viewpoint that mass must move along edges and the transportation should be accounted in terms of edge flows, which are really just signed mass values on each oriented edge. In this formulation, the penalty lies on |J(e)|^p, so transportation flows can often be non-sparse, or occur along weighted combinations of shortest and non-shortest paths. Comparisons of optimal flows for values p=1, 2 in an example setting are shown in <ref>.

	A slightly different way to formulate <ref> would be through E” with a nonnegativity constraint. That is, let
			
    𝒥”(α,β)    = {J∈ℓ_2(E”) : J≥ 0, BJ=α-β,  }.

		and then,
			
    ℬ_p(α,β)    = inf{J_w,p : J∈𝒥”(α,β)}.

		It is not hard to show that these two formulations arrive at the same value. <ref> simply allows for signed flows; for example, a flow value of -1 along the oriented edge (1, 2) represents one unit of mass transported from node 2 to node 1. We opt to drop the nonnegativity constraint for this paper and focus on the former formulation as it simplifies the calculations slightly, but we may have occassion to use this formulation instead.
	

	In the Proposition below, we summarize two relevant properties of 𝒲_p and ℬ_p.

	
		Let α,β∈𝒫(V).
		
			
  * The objective functions for both 𝒲_p(α,β) and ℬ_p(α,β) are convex for all 1≤ p<∞, and the infima are always achieved.
            
  * The objective for ℬ_p(α,β) is strictly convex for p>1. The infimum is always achieved and is unique in the case p>1.
			
  * 𝒲_p and ℬ_p(α,β) define metrics on 𝒫(V). 
		
	

    Proofs of items (i) and (ii) are straightforward since the objectives are defined either by affine functions or common norms whose convexity properties are well-understood. For (iii), a proof for 𝒲_p can be found in <cit.>, and the case of ℬ_p is relatively straightforward.
    
	

§ GENERAL PROPERTIES OF THE BECKMANN PROBLEM

	
	In this section, we look at a few interesting properties of the Beckman problem ℬ_p(α,β) for general 1≤ p<∞. In <ref> we derive the Lagrangian dual to the p-Beckmann problem, in <ref> we obtain bounds that relate certain Beckmann and Wasserstein distances, and finally, in <ref>, we work through some example results on paths and trees to highlight examples of Wasserstein and Beckmann distances.
	
	

 §.§ Duality theory for the p-Beckmann problem

	
	In this subsection we study the convex duality of the p-Beckmann problem. In <cit.> Alamgir and Luxburg introduced formulations of p-resistances, that is, generalizations of effective resistance obtained by Beckmann-type min cost flow problems penalized by a p-norm. They studied the primal and dual forms of these problems, and thereby obtained so-called potential-based formulations of the resistance problems. As described in <ref>, the results presented in this section may be considered extensions of the results concerning p-resistances to the case of general probability measures; namely, we understand the p-resistance between generic nodes i, j to be p-Beckmann distance between Dirac measures at i, j. In <ref> we obtain the dual formulation of the p-Beckmann optimization problem, and thus generalize the Kantorovich-Rubenstein duality on graphs to all 1≤ p<∞ (see <ref>).
    
    
	
	
		Let α,β∈𝒫(V) and 1 < p <∞ be fixed. Let q be the conjugate of p so that 1/p + 1/q= 1. Then the dual to the p-Beckmann problem, given by
			
    ℬ_p(α,β)    = inf{J_w,p : J∈ℓ(E'), BJ = α-β},

		is given by the maximization problem
			
    ℬ_p(α,β)    =sup{ϕ^T(α-β) : ϕ∈ℓ(V), B^Tϕ_w^1-q, q≤ 1 },

		and strong duality holds. In the special case of p=1 and q=∞, we have
			
    ℬ_1(α,β)    = inf{J_w,1 : J∈ℓ(E'), BJ = α-β} = sup{ϕ^T(α-β) : ϕ∈ℓ(V), B^Tϕ_w^-1, ∞≤ 1 }

	

	We defer the proof to <ref>.















	If we specialize to the case where the edge weights are all unit value; i.e., G is unweighted, then the weights can be dropped from the preceding results. In particular, the dual norm to ·_p is simply ·_q and 
			
    ℬ_p(α,β) = inf{J_p : J∈ℓ(E'), BJ = α-β} = sup{ϕ^T(α-β) : ϕ∈ℓ(V), B^Tϕ_q≤ 1 }

	
	
	The case p=1 and q=∞ is the conventional Kantorovich-Rubenstein duality on graphs <cit.>, which can be more conventionally written as
		
    inf{∑_e∈ E' |J(e)|w_e : J∈ℓ(E'), BJ = α-β} = sup{ϕ^T(α-β) : ϕ∈ℓ(V), ϕ_Lip≤ 1 }

		where the Lipschitz seminorm ϕ_Lip = B^Tϕ_w^-1,∞ on ℓ(V) is defined by
        
		
    ϕ_Lip = max|ϕ(i) - ϕ(j)| /d(i, j) = max_i∼ j|ϕ(i) - ϕ(j)|/w_ij.

        Note also that this duality occurs in a homophilous manner for Beckmann's minimal flow problem in the continuous setting <cit.>, which is the namesake for the minimal flow formulations on graphs.
        
	
	
	
	

 §.§ Comparing p-Beckmann and p-Wasserstein distances

	
	In this section, we establish some bounds between ℬ_p and 𝒲_p in certain situations. <ref> introduces a technique for obtaining an edge flow from a given coupling with some control on the cost, and uses it to derive elementary bounds for general p; and in <ref> we focus on the case of 𝒲_1 and obtain sharp estimates related to ℬ_2, the latter of which is of particular interest in the ensuing sections.
	
	

  §.§.§ From coupling to flow

	
	An essential ingredient in comparing the optimal values of transportation metrics is having the ability to convert a coupling π to a feasible edge flow J, and vice-versa. In this subsection we will make use of the bi-oriented edge flow formulation of ℬ_p, which was discussed in <ref>.
	
	Let i, j∈ V. Let P∈𝒫(i, j) be a path between i, j. Write P = (i = i_0, i_1,…, i_k = j). Define 
			
    I_P = ∑_ℓ=0^k-1δ_(i_ℓ , i_ℓ+1)∈ℓ(E”).

	
	
	
		Let α,β∈𝒫(V), 1≤ p <∞, and suppose π∈Π(α,β). For each pair of nodes i, j∈ V with i≠ j, let P_ij∈𝒫(i, j) be a fixed path between i, j. Define the edge flow J^π = J^π, (P_ij)∈ℓ(E”) by the expression
			
    J^π = ∑_i≠ j∈ Vπ_ij I_P_ij
 
		Then J^π is a feasible flow for ℬ_p(α,β), i.e., J^π≥ 0 and BJ^π=α-β, and the following estimate holds:
			
    J^π_w,p≤∑_i, j∈ Vπ_ijI_P_ij_w,p.

	
	
	We include a proof in <ref>. We can use <ref> to derive at least two estimates which relate ℬ_p to 𝒲_k, p for different choices of k. The first is somewhat brutal, but is informative for p small.
	
	
		Let α,β∈𝒫(V), 1≤ p < ∞, and assume that w_ij≥ 1 for each i, j∈ V such that i∼ j. Then
			
    ℬ_p(α,β)≤𝒲_p(α,β)^p

	
	
	
		Let π∈Π(α,β) be an optimal coupling for 𝒲_p(α,β). For each i, j, let P_ij∈𝒫(i, j) be a choice of path which is minimal in the sense of d_1, and let J^π = J^π, (P_ij) be the feasible flow as in <ref>. Then
			
    ℬ_p(α,β) ≤J^π_w,p≤∑_i, j∈ Vπ_ijI_P_ij_w, p≤∑_i, j∈ Vπ_ij d_1(i, j)^p =𝒲_p(α,β)^p

		since ·_w, p≤·_w, 1^p provided w_ij≥ 1.
	
	
	If we allow for some flexibility when picking the underlying metric k of 𝒲_k, p, then we can provide a tighter estimate when 𝒲 is large.
	
	
		Let α,β∈𝒫(V), 1< p < ∞ and let q be the conjugate of p such that 1/p + 1/q = 1. Then
			
    ℬ_p(α,β)≤ n^2/q𝒲_d_p, p(α,β).

	
	
	The proof of <ref> follows mainly from Hölder's inequality and <ref>; we provide the proof in <ref>. Note in particular that when p=2, we have q=2 and thus it holds ℬ_2(α,β)≤ n 𝒲_d_2, 2(α,β).
    

	

  §.§.§ Estimating from the p=1 case


	When p=1 and the choice of metric is k=d_1 there is a direct overlap between ℬ_1 and 𝒲_1 <cit.>.
	
	
		Let α,β∈𝒫(V). Then 
			
    ℬ_1(α,β) = 𝒲_1(α,β).

	
	
	For a proof of this result, we recommend the exposition in <cit.>. Note that for p≥ 1 and J∈ℓ(E'), we have that since ·_p≤·_1 in general,
		
    J_w, p =  W^1/pJ _p ≤ W^1/pJ _1 ≤J_w, 1(max_e∈ E' w_e^1/p-1).

	Thus if C_w,p = max_e∈ E' w_e^1/p-1 we have J_w, p≤ C_w,pJ_w, 1. If J^∗ is an optimal flow for ℬ_1(α,β), we have
		
    ℬ_p(α,β)≤J^∗_w, p≤  C_w,pJ^∗_w, 1 = C_w,pℬ_1(α,β)

	Therefore, we have the following two estimates for any choice of α,β and p≥ 1.
	
	
		Let α,β∈𝒫(V) and 1≤ p <∞. Then the following estimates hold:
			
    ℬ_p(α,β)   ≤ C_w,pℬ_1(α,β)
    ℬ_p(α,β)   ≤ C_w,p𝒲_1(α,β)

	
	
	Moreover, via similar logic as before, since ·_1 ≤ m^1 - 1/p·_p on ℓ(E'), where we recall m = |E|, we have for J∈ℓ(E'),
		
    J_w, 1 = W J_1 ≤ m^1 - 1/p W J_p ≤ m^1 - 1/pW1_E'_∞^p-1/pJ_w,p

	Therefore if we set C_w,m,p = m^1 - 1/pW1_E'_∞^p-1/p, we have J_w, 1≤ C_w,m,pJ_w,p and using a similar argument as before, the following corollary holds.
	
	
	Let α,β∈𝒫(V) and 1≤ p <∞. Then the following estimates hold:
		
    ℬ_1(α,β)   ≤ C_w, m, pℬ_p(α,β)
    𝒲_1(α,β)    ≤ C_w, m, pℬ_p(α,β)

	
	
	An application of this allows us to relate 2-Beckmann distance to 𝒲_1, which we state as a corollary below.
	
	
		Let α,β∈𝒫(V). Then 
			
    ℬ_2(α,β) ≤ C_w,2𝒲_1(α,β) ≤ C_w, 2C_w, m, 2ℬ_2(α,β).

		If in particular the graph is unweighted, then we have that C_w, 2 = 1 and C_w, m, 2 = m^1/2, so that
			
    ℬ_2(α,β) ≤𝒲_1(α,β) ≤ m^1/2ℬ_2(α,β).

	

    Note that the bound in <ref> is sharp without making additional assumptions on the structure of G or α,β. Suppose G is a path on n vertices with m=n-1, then the upper and lower bounds are achieved, respectively, when we have:
       		
       			
  * α=δ_1, β=δ_n, so that ℬ_2(δ_1, δ_n) = √(n-1) and 𝒲_1(δ_1, δ_n) = n-1 so 𝒲_1 = m^1/2ℬ_2.
       			
  * α=δ_1 and β=δ_2, so that ℬ_2 = 𝒲_1 = 1.
       		
        Proofs of these statements could be obtained directly or by way of <ref> and <ref>, which we discuss in the next subsection. This leads to an as-yet open question: If one constrains the structure of G or α,β, can the bound in <ref> be improved?
    

	

 §.§ Worked examples: Paths and trees

 
    
 
	In this subsection, we specialize to the case of the path graph P_n and derive explicit formulas for 𝒲_p and ℬ_p so as to illustrate their connections to optimal transportation in the setting of the real line and to elucidate their differences when p > 1. We also take a look at ℬ_p on tree graphs, and obtain a closed-form solution for all 1≤ p<∞ which generalizes the known formula for 1-Wasserstein distance on trees. Proofs of <ref>, <ref>, and <ref> are all included in <ref>.
	
	Let P_n be the unweighted path graph with n nodes ordered from 1 to n, and oriented edges E' of the form (i, i+1) for 1≤ i ≤ n-1. For α∈𝒫(V(P_n)), let α = ∑_i=1^n α(i) δ_i be the empirical measure on ℝ induced by α where, with a slight abuse of notation, δ_i is the Dirac measure on ℝ at x=i.
    
    Let F_α be the cumulative distribution function of α, 
    
    and let F_α^-1 be its psuedo-inverse or quantile function, i.e., F_α^-1(t) = inf{x : F_α(x) ≥ t} for t∈ [0, 1]. The proposition below is an adaptation of the classical inverse cdf result for 𝒲_p on ℝ to the case of path graphs; for more detail, see, e.g., <cit.>.
	
	
		Let α,β∈𝒫(V(P_n)) and 1≤ p <∞. Then the p-Wasserstein distance 𝒲_p(α,β) on P_n is given by
			
    𝒲_p(α,β)^p = ∫_0^1 | F_α^-1(t) - F_β^-1(t)|^p  t.

	
	
	To characterize p-Beckmann distance on P_n, we introduce a bit more notation. For α∈ℙ(V) and an oriented edge e=(i, i+1)∈ E', let K_α∈ℓ(E') be the “edge cdf," i.e.,
		
    K_α(i, i+1) = ∑_j ≤ iα(j).

    Note that, naturally, K_α(i, i+1) = F_α(i).
    
	
	
		Let α,β∈𝒫(V(P_n)) and 1≤ p <∞. Then the p-Beckmann distance ℬ_p(α,β) on P_n is given by
			
    ℬ_p(α,β)^p = K_α - K_β_p^p = ∫_ℝ | F_α(t) - F_β(t) |^p  t.

	
	
	<ref> and <ref> provide an alternate, albeit circuitous, proof of <ref> on paths by applying a change of variables. This also demonstrates why we do not expect the two families of metrics to overlap for p>1. A worked example of <ref> and <ref> in action is illustrated in <ref>. 
 
    It is also intriguing to note that much in the same way that 𝒲_1, and hence ℬ_1, is determined explicitly on trees (a result which itself has been leveraged in other applications such as computer vision and natural language processing <cit.>), ℬ_p can be determined in an similar manner for general p on trees.
	
	
		Let T=(V, E, w) be a weighted tree, α,β∈𝒫(V), and fix 1≤ p <∞. For an oriented edge e = (i, j)∈ E', define a generalized version of <ref> by
			
    K_α(e= (i, j)) = ∑_k ∈ V^∗(i; e)α(k),

		where V^∗(i; e)⊂ V is the set of nodes belonging to the subtree with root i obtained from T by removing the edge e (and similarly for K_β). Then it holds
			
    ℬ_p(α,β) = K_α-K_β_w,p.
 
	
 
	

§ THE 2-BECKMANN PROBLEM: THREE PERSPECTIVES

	
	This section is focused entirely on ℬ_2. The foundational premise of this section is that ℬ_2, while ostensibly a simple least squares optimization problem, actually posesses a rich supply of connections and contexts to other notions that already exist in graph theory. In <ref> we view ℬ_2 as an effective resistance metric between measures, and connect it to the simple random walk on G and optimal stopping times for measures. In <ref>, we view ℬ_2 as a negative Sobolev norm for functions on graphs, and thereby obtain a Benamou-Brenier type formula for ℬ_2. Lastly, in <ref>, we view ℬ_2 as a linearized variant of 𝒲_2, prove a convex separation result, and showcase its potential application in unsupervised learning problems for graph data.
	
	

 §.§ Effective resistance between measures

	
	In this subsection we consider ℬ_2 from the perspective of effective resistance. We begin with some supplemental background on effective resistance which will be needed, and then proceed with discussion of the main results.
	
	

  §.§.§ Background on effective resistance

	
	Let i, j∈ V be any two nodes. The effective resistance between i, j, denoted r_ij, is given by the formula
			
    r_ij = (δ_i - δ_j)^T L^† (δ_i - δ_j).

	
	
	The effective resistance between the nodes of a graph has many different properties and representations, many of which intersect with the simple random walk on G. We give its definition and some of its properties below.
	
	The simple random walk on G is the Markov chain (X_t)_t≥ 0 on the state space of nodes V with transition probability matrix D^-1A; that is,
			
    ℙ[X_t+1 = j | X_t = i] = w_ij/d_i    if i∼ j
    
    				0   otherwise .

	
	
	Recall that X_t will admit a stationary distribution whenever G is connected and will be ergodic whenever G is also non-bipartite. The  stationary distribution, which we denote ρ, is always proportional to the degree at each node. A useful tool is the vertex hitting time
		
    T^s_j = inf{t≥ s: X_t = j},s≥ 0, j∈ V.

    Also, let (G) = ∑_i∈ V d_i denote the volume of the graph G.
	 
	
		Let i, j∈ V be any two nodes. Then the following observations hold
			
				
  * (i, j)↦ r_ij is a metric on V,
				
  * (i, j)↦ r_ij^1/2 is a metric on V,
				
  * r_ij = inf{ f^T L f : Lf = δ_i - δ_j },
				
  * r_ij = L^-1/2 (u_i - u_j)_2^2, where u_i is the i-th orthonormal eigenvector of the graph Laplacian L.
				
  * r_ij = 1/vol(G)(𝔼[T^0_i : X_0 = j] + 𝔼[T^0_j : X_0 = i]), i.e., r_ij is proportional to the commute time between i,j.
			
	
	
	For proofs of these results, see, e.g., <cit.>.
	
	

  §.§.§ 2-Beckmann as measure effective resistance

	
	We begin with the observation that ℬ_2(α,β)^2 can be realized as a type of effective resistance for measures. We make this precise in the theorem below.
	
	
		Let α,β∈𝒫(V). Then 
			
    ℬ_2(α,β)^2 = (α-β)^T L^†(α-β).

	

	
		We can write
			
    ℬ_2(α,β)^2    = inf{J_w, 2^2 : BJ=α-β}
       = inf{ W^1/2 J_2^2 : BJ=α-β}
       = inf{ S_2^2 : (BW^-1/2)S=α-β}

		with S = W^1/2J. Recall that for matrices X, x it holds that X^† x=argmin{y_2 :  Xy =x} when x belongs to the column space of X <cit.>, and thus we have
			
    ℬ_2(α,β)^2    =  (BW^-1/2)^† (α-β)_2^2
       = (α-β)^T ((BW^-1/2)^†)^T(BW^-1/2)^† (α-β)
       = (α-β)^T L^† (α-β)

		since L = BWB^T. 
	
	
	One immediate consequence of <ref> is a Brenier-type theorem for graphs, which we state as a corollary below. Based on the discussion following <ref>, in this theorem we opt to use the notation B^T = ∇ to achieve homophily with the continuous setting.
	
	
		Let α,β∈𝒫(V). The unique optimal edge flow J^∗ for  ℬ_2(α,β) can be written J^∗ = ∇ f where f is a mean-zero potential function of minimal ℓ_2 norm satisfying the Poisson-type equation Lf = α-β. The optimal cost is ∇ f_w, 2. 
	
	
	
		It is enough to write (α-β)^T L^† (α-β) = (α-β)^T L^† L L^† (α-β) and L = BWB^T. The claim then follows from the proof of <ref>.
	
	
	Based on the <ref>, we introduce the following definition and notation of effective resistance for measures.
	
	Let α,β∈𝒫(V). Then we define the measure effective resistance by 
			
    r_αβ = (α-β)^T L^†(α-β).

	

    In the subsequent sections, we will explore the properties of r_αβ in its own right, with particular focus on its relationship to optimal stopping rules of the simple random walk.
    
    
	

  §.§.§ Measure effective resistance and optimal stopping rules

	
	We have already argued in <ref> that measure effective resistance can be understood as a type of flow-based optimal transportation distance between the measures. The goal of this section is to essentially set that aside and think about r_αβ in its own right. A natural approach along these lines is to take <ref> as a starting point. Immediately, r_αβ has several analagous properties, which we state below without proof.
	
	
		Let α,β∈𝒫(V). Then the following statements hold:
			
				
  * (α,β)↦ r_αβ^1/2 is a metric on 𝒫(V),
				
  * r_αβ = inf{ f^T L f : Lf = α-β},
				
  * r_αβ = L^-1/2(α-β)_2^2.
			
	
	
	To go deeper, we pose the following question along the lines of <ref>(v): To what extent can r_αβ be understood as a generalized commute time between the measures α,β? To answer this question, we need to use results that concern stopping rules for the simple random walk.
 
	A stopping rule is a map Γ that associates to each finite path ω = (X_0, X_1, …, X_k) on G a number Γ(ω) in [0, 1]. We can think of Γ(ω) as the probability that we continue a random walk given that ω is the walk so far observed. Alternatively, Γ can be considered a random variable taking values in {0, 1, 2, …} whose distribution depends only on the steps (X_0, X_1, …, X_Γ).
	
	
	The mean length 𝔼[Γ] is the expected duration of the walk. If 𝔼[Γ] < ∞ then the walk stops almost surely in finite time, so we define X_Γ to be the position of the random walk at the stopping time. Having defined stopping rules, we can define the generalized hitting time between α, β.
	
	Let α,β∈𝒫(V). The access time H(α,β) is defined as
			
    H(α,β) = inf{𝔼[Γ | X_0∼α] : 𝔼[Γ] < ∞, and  X_Γ∼β}.

		where for any random variable Y on V, we say Y ∼α if ℙ[Y = i] = α(i) for i∈ V. In other words, H(α,β) is the minimum mean length of walks that originate with distribution α and terminate according to a stopping rule that achieves distribution β at stopping time. If Γ achieves the inf in H(α,β), then Γ is said to be an optimal stopping rule.
	
	
	It is not hard to see that the set of feasible stopping rules in <ref> is nonempty. The so-called “naïve" stopping rule Γ_n can be obtained from the following construction: at the beginning of the random walk, sample j∼β, and stop the walk when X_Γ_n = j. It is readily verified that X_Γ_n∼β, and that
		
    𝔼[Γ_n] = ∑_i, j∈ Vα_iβ_j H(i, j)

		where for i, j∈ V, the hitting time H(i, j) is defined by H(i, j) = H(δ_i, δ_j) (or, the mean number of steps to reach j from i).
	

    
	
	Several examples of optimal stopping rules for general distributions are given in <cit.>, and in particular, it is shown thereby that any such α,β admit an optimal stopping rule provided G is connected. 
	
	The access time H(α,β) has several notable properties, which we summarize in the Proposition below. The proofs of these results can be found in <cit.>.
	
	
		Let α,β∈𝒫(V). Then the following facts hold for H(α,β):
			
				
  * When β is concentrated at a node j∈ V, it holds H(α,j) = ∑_i∈ Vα_i H(i, j),
				
  * H(α,β) = max_j∑_i∈ V (α_i-β_i) H(i,j),
				
  * H(α,β) is convex in both of its arguments; namely if α,α', β,β'∈𝒫(V) and c∈ [0, 1], we have
					
    H(c α + (1-c)α', β)    ≤ cH(α, β) + (1-c)H(α', β)
    
    						H(α, c β + (1-c)β' )   ≤ cH(α,β) + (1-c)H(α, β')

			
	
	
	Beveridge <cit.> also established a connection between the entries of the psuedoinverse of the graph Laplacian (treated as a discrete Green's function) and the access time between nodes and the stationary distribution of the random walk. We reproduce that formula here as a proposition, in a modified form to match our conventional choice of Laplacian matrix.
	
	
		Let ρ be the stationary distribution of the simple random walk on V. Then the i, j-th entry of L^† can be expressed as 
			
    (L^†)_ij = 1/(G)( H(ρ, j) - H(i, j)),

		where (G) = ∑_i∈ V d_i.
	
	
	For α∈𝒫(V), define H_α∈ℓ(V) by H_α(i) = H(α, i). We are ready to state our main result connecting r_αβ to the access time.
	
	
		Let α,β∈𝒫(V). Then 
			
    r_αβ = -1/(G) (α-β)^T (H_α - H_β).

		Or, in alternative expanded forms,
			
    r_αβ   = -1/(G)∑_i∈ V(α_i - β_i) (H(α, i) - H(β, i))
       = -1/(G)∑_i, k∈ V(α_i - β_i)   (α_k - β_k) H(i, k).

	
	
	We provide the proof in <ref>. The reason we term <ref> a generalized commute time formula is that in the case where α,β are Dirac measures, this recovers the result that r_ij∝ (H(i, j) + H(j, i)), since we have
        
    r_ij    = -1/(G) (δ_i - δ_j)^T(H_δ_i - H_δ_j)
       = -1/(G) (H(i, i) - H(i, j) + H(j, j) - H(j, i)
       = 1/(G)(H(i, j) + H(j, i)).

    <ref> shows, however, that in the more general case of measures r_αβ is not exactly the same as a measure commute time H(α,β) + H(β,α). <ref> can be used to establish a relationship between the measure commute time and the measure resistance, as we show below.
	
	
		Let α,β∈𝒫(V). Then r_αβ satisfies the following two inequalities:
			
    r_αβ   ≤2/(G)max{H(α,β), H(β,α)}
    
    				r_αβ   ≤1/(G) (H_n(α,β) + H_n(β, α))

		where H_n(α,β) = 𝔼[Γ_n] (resp. H_n(β,α)) is the expected duration of the naïve stopping rule described in <ref> with initial distribution α (resp. β) and stopping node sampled from β (resp. α).
	
	
	The proof of <ref> is located in <ref>. Note that the second half of the proof of <ref> actually contains a related observation regarding ℬ_2 and 𝒲_r, 1, which we state as another corollary below, without proof.
	
	
		Let α,β∈𝒫(V), and let r be the effective resistance metric on V. Then we have
			
    ℬ_2(α,β) ≤𝒲_r, 1(α,β)^1/2.

	 
	
	A final observation that one can glean from the random walk and access time setting is an interpretation of the transportation potential f = L^†(α-β) mentioned in <ref>.
	
	
		Let α,β∈𝒫(V), and let Γ be any stopping rule such that 𝔼[Γ] < ∞ and X_Γ∼β. Let f_Γ be the degree-weighted exit frequencies
			
    f_Γ(i) = 𝔼[∑_t=0^Γ - 11/d_iδ_{X_t = i}(i) : X_0∼α].

		Then L f_Γ = α-β. In particular, the transportation potential L^†(α-β) is given by the mean normalized function f_Γ - 1/n f_Γ^T 1_n for any suitable stopping rule Γ.
	
	
	This proposition follows from applying the so-called conservation equation for exit frequencies, mentioned in, e.g., <cit.>. Another way of putting it is that the unique optimal flow for ℬ_2(α,β) can be realized as the gradient of the degree-normalized exit frequencies of any finite-mean stopping rule which achieves a stopping distribution of β, having been intialized at α (since, moreover, as discussed in <cit.>, any two degree-normalized exit frequencies for different stopping times differ by a constant, and thus have the same gradient). 

	

 §.§ Sobolev norms and a graph Benamou-Brenier formula

	
	The relationship between Sobolev spaces and optimal transportation have been explored at length in the continuous setting <cit.>. The goal of this subsection is to present the 2-Beckmann problem as an optimization problem on a negative Sobolev space for functions defined on graphs. In this section we follow the exposition laid out in <cit.> in the continuous setting by analogy to the discrete setting. As it turns out, there are at least two results that achieve strong homophily between the graph and continuous setting.
	
	

  §.§.§ Background from the continuous setting

	
	Recall that if f:ℝ^n→ℝ is a function with a square integrable derivative ∇ f in the weak sense and μ is a Borel probability measure which is absolutely continuous with respect to the Lebesgue measure x; i.e. so that μ = g  x for a density function g, we can define the Sobolev-type seminorm ·_Ḣ^1(μ)^2 by
		
    f_Ḣ^1(μ)^2 = ∫_ℝ^n∇ f_2^2 μ.

	The dot Ḣ^1(μ) serves to distinguish ·_Ḣ^1(μ)^2 from a true Sobolev norm, which include a contribution from ·_L^2. We can then define the possibly infinite dual norm to f_Ḣ^1(μ)^2, denoted ·_Ḣ^-1(μ) by the following, for any x- absolutely continuous signed measure ν = h  x:
		
    h x _Ḣ^-1(μ) = sup{∫_ℝ^n f h dμ : f_Ḣ^1(μ)≤ 1 }.

	This setup leads to two results (among others) in the continuous setting, which are of interest in the graph setting as well. 
	
	
		If μ is a Borel probability measure on ℝ^n and μ is an infintesimally small perturbation of μ, then
			
    𝒲_2(μ,μ + μ ) = μ_Ḣ^-1(μ) + o(μ).

	
	
	The second result which interests us is the Sobolev form of the Benamou-Brenier formula, which we state below.
	
	
		Let μ,ν be Borel probability measures on ℝ^n. Then it holds:
			
    𝒲_2(μ,ν) = inf{∫_0^1 μ_t_Ḣ^-1(μ_t) : μ_0 = μ, μ_1 = ν}.

	
	
	

  §.§.§ Graph sobolev norms

	
	Next we transition to defining the appropriate discrete analogues of ·_Ḣ^1(μ) and ·_Ḣ^-1(μ) on the graph G.
	
	Let f, g∈ℓ(V). We define the graph Sobolev seminorm ·_Ḣ^1(V) by the equation
			
    f_Ḣ^1(V)^2 = ∑_(i, j)∈ E' w_ij |∇ f(i, j)|_2^2 = ∑_(i, j)∈ E'w_ij |f(i) - f(j)|_2^2.

		We define the (possibly infinite) dual graph Sobolev norm ·_Ḣ^-1(V) by the supremum
			
    g_Ḣ^-1(V)^2 = sup{ f^T g: f_Ḣ^1(V)≤ 1 }.

	
	
	It is useful to note that for mean zero functions, ·_Ḣ^1(V) and ·_Ḣ^-1(V) will be true norms (in particular, the former will be definite and the latter will be finite). The first proposition in the section elucidates the relationship between graph Sobolev norms and the graph Laplacian matrix.
	
	
		Let f, g∈ℓ(V). Then the following hold:
			
				
  * f_Ḣ^1(V)^2 = f^T Lf.
				
  * If 1^T g = 0, then g_Ḣ^-1(V)^2 = g^T L^† g.
			
	 
	
	
		The first claim is straightforward. The second requires a closer look. We can work backwards and use the proof of <ref> to obtain
			
    (g^T L^† g)^1/2   =  L^-1/2 g _2 
       =  inf{ f_2: L^1/2 f = g }
       = sup{ f^T g:  L^1/2 f_2 ≤ 1 },

		where L^{-1/2, 1/2 } = U Λ^{-1/2, 1/2} U^T as in <ref>. But L^1/2 f_2^2 = f^T L f = f_Ḣ^1(V)^2, so the claim follows.
	
	
	

  §.§.§ 2-Beckmann as a negative Sobolev distance

	
	Since our definitions have been set up and explained properly, we can derive the following result, which follows directly from <ref> and <ref>. We consider this a graph analogue of <ref> for 2-Beckmann distance. 
	
	
		Let α,β∈𝒫(V). Then
			
    ℬ_2(α,β) =α - β_Ḣ^-1(V).

		Thus in particular if β = α + α for a small perturbation α, then ℬ_2(α,α + α) = α_Ḣ^-1(V).
	

	Note that since the Ḣ^-1 norm depends on the operator L^-1/2, the Beckmann distance ℬ_2(α,α + α) can be articulated in terms of the spectral coefficients of α; to wit, if α = ∑_ℓ=2^n c_ℓ u_ℓ is the spectral decomposition of a mean-zero perturbation α∈ℓ(V) in terms of the eigenvectors u_ℓ of L, then 
		
    ℬ_2(α,α + α)^2 = L^-1/2α_2^2 =∑_ℓ=2^n c_ℓ^2/λ_ℓ.

	Therefore in general, one has ℬ_2(α,α + α)≤λ_2^-1/2α_2; but if the spectral properties of α are known, this estimate may be sharpened.









	
	The next result which leverages the Sobolev norm perspective on graph optimal transport can be considered a graph analogue of the Benamou-Brenier formula. Let μ_t∈ℓ(V) for each t∈ [0, 1].
	We say μ_t∈ C^1([0, 1]) if the map t↦μ_t : [0, 1]→ℓ(V) is continuously differentiable as a map from [0, 1] to ℝ^n. Moreover, when μ_t admits a derivative, we write μ_t = . / sμ_s|_s=t.
	
	
		Let α,β∈𝒫(V). Then we have
			
    ℬ_2(α,β)^2 = inf{∫_0^1 μ_t _Ḣ^-1(V)^2  t: μ_t∈ C^1([0, 1]), μ_0 = α, μ_1 = β}.

	
	
	See <ref> for a proof, which is based on the proof appearing in the lecture notes <cit.>. Note that the infimum is achieved by a linear line segment connecting α to β, which underlines the idea that ℬ_2 should be treated as a linearized variant of Wasserstein distance. We explore the linearization angle and its implications on graph learning tasks in the next section.
	
	

 §.§ A linearized optimal transportation distance

	
	A third and final perspective on 2-Beckmann distance is from the world of clustering and classification of graph data. We begin with a background discussion from the continuous setting, and then explore convex separation properties for graph data as well as an application involving a handwritten digit dataset.
	
	

  §.§.§ Background from the continuous setting

	
	A typical classification scenario usually consists of some data {x_i}⊂ℝ^n which one wishes to separate into classes or clusters either without prior knowledge of the data labels (unsupervised) or with some information about the class labels of the dataset (semi-supervised or supervised). In many applications, e.g., <cit.>, the data x_i can often occur not as vectors in ℝ^n but as distributions μ_i on ℝ^n.
	
	In such scenarios, Wasserstein distance, formulated in the Monge sense as 
		
    𝒲_2(μ,ν)^2 = inf_T:T_♯μ = ν∫_ℝ^nT(x) - x_2^2 dμ(x),

	where T_♯μ is the push-forward of μ under the map T and the inf is over all such T exhibiting desired regularity, arises as a natural distance metric between probabilitiy measures on ℝ^n. However, this metric can be slow to compute at scale, and in the discrete setting can even suffer from being undefined. 
	
	One approach to resolving this issue involves the usage and study of Linearized Optimal Transport embeddings <cit.>, which can be summarized as follows. Given data measures {μ_i}, fix a reference measure μ, and define T_i = argmin_T:T_♯μ = μ_i∫_ℝ^nT(x) - x_2^2 dμ(x). Then the measures are featurized in the form of the maps T_i:ℝ^n→ℝ^n and 2-Wasserstein distance between the data measures is then approximated by T_i - T_j_L^2.

	Theoretically, the ability of Linearized Optimal Transport in approximating 2-Wasserstein distance can be captured by the following two results, the informal statements of which we restate from <cit.>.
	
	
		If 𝒫 = {μ_i} are ϵ-perturbations of shifts and scalings of μ, and 𝒬 = {ν_i} are ϵ-perturbations of shifts and scalings of ν, and 𝒫 and 𝒬 have small minimal distance depending on ϵ and satisfy a few technical assumptions, 𝒫 and 𝒬 are linearly separable in LOT embedding space.
	

	
		If μ,ν are ϵ-perturbations of shifts and scalings of one another, then
			
    𝒲_2(μ, ν) ≤𝒲_2^LOT (μ, ν) ≤𝒲_2(μ, ν) + C_σϵ + C_σ'ϵ^1/2.

		In particular, when ϵ=0 and μ,ν are shifts and scalings of each other, LOT is an isometry.
	

	It is also worth noting that in the situation where μ_i is some shift or scaling of μ, i.e., μ_i = T_♯μ for an affine combination of scaling and translation T on ℝ^n, then the LOT embedding of μ_i will be the map T itself; i.e., T is optimal in the Wasserstein cost.
	
	

  §.§.§ Linear separtation for 2-Beckmann distance

	
	To set up the problem, suppose we have data sampled from C classes, and each class is associated to a canonical distribution α_i ∈𝒫(V) so that our data may be represented in the form 𝒟 = {T α_i}_i=1,…, C; T∈𝒯. That is, each data point is of the form T α_i where α_i is the underlying distribution of the class and the map T:𝒫(V)→𝒫(V) represents a perturbation of the class distribution among valid perturbations 𝒯 (to be determined).
	
	For example, consider a hypothetical dataset of images with resolution k×ℓ. In this setting, G is the k×ℓ lattice graph, and after normalization each image can be understood as a distribution on G. In our data model, we assume that each image may be expressed as a perturbation of a canonical distribution associated to the corresponding image class.
	
	In the unsupervised setting, we simply have 𝒟 and we wish to use some clustering technique, e.g., k-means or principal geodesic analysis, associated to a metric on 𝒫(V). In the supervised setting, we have data equipped with labels 𝒟' = {(Tα_i, y_i)}_i=1,…, C; T∈𝒯 so that each label depends only on the class distribution and we seek to build some classifier, ideally linear, which facilitates predictive analysis.
	
	We treat ℬ_2 as a linearization of 𝒲_2 on graphs. In this setting, the featurization of each canonical class distribution is given by α_i↦ L^-1/2α_i and similarly for the individual samples Tα_i. 
	
	For a set of measures 𝒜 = {α_1,…,α_k}⊂𝒫(V), define the mutual support by
        
    ℳ𝒮(𝒜) = ⋂_i=1^k supp(α_i)⊂ V,

    where supp(α_i) is the subset of nodes in V on which α > 0.
	
	
		Let α_1≠α_2∈𝒫(V) be the canonical class distributions for a binary distributional dataset 𝒟 = {T α_i}_i=1,2; T∈𝒯. Assume ℳ𝒮({α_1,α_2})≠∅ and let 𝒯 be the set of affine perturbation maps defined by
			
    Tα = α +  T

		where T∈ℓ(V) is vector which satisfies the three conditions:
			
				
  * T is mean zero, or 1_n^T T = 0;
				
  * ( T) ⊂ℳ𝒮({α_1,α_2}); and
				
  * T_1 < δ, where
					
    δ = min{1/3α_1 - α_2_2,  min_i=1, 2
    k∈ℳ𝒮({α_1,α_2}) α_i(k)} > 0.

			
		Then the sets 𝒜_1 = {Tα_1}_T∈𝒯 and 𝒜_2 = {Tα_2}_T∈𝒯, which belong to the metric space (𝒫(V), ℬ_2), can be isometrically embedded into ℓ_2(V) under L^-1/2 so that their images are linearly separable as convex sets in ℝ^n.
	
	
	See <ref> for the proof. Note that it is not hard to extend <ref> to the case of C>2 classes, by simply requiring that the perturbations be smaller than the mutual distances between the canonical class distributions, up to an appropriate constant, so that in the data model there are no overlapping perturbations. The main point of this result is that L^-1/2 as an embedding of measures, is well-behaved; in the sense that as long as the original collections of measures were separable, L^-1/2 will not imbue any pathology.	

	

  §.§.§ Experimental results


    
 
    
	
    
	
	To further explore the implications of our results on learning with graph data, we conducted two experiments which illuminate some of the potential applications and underlying properties of 2-Beckmann distance between measures. Both examples make use of the “Pen-Based Recognition of Handwritten Digits" dataset <cit.>, accessed and processed via the Sklearn Python package <cit.>. The dataset contains 1797 grayscale handwritten digits of resolution 8× 8 and pixel values ranging from 0-15. The reason for choosing this dataset is that the 2-Wasserstein kernel is computationally prohibitive at scale without approximation, so rather than regularize the Wasserstein metric, we opt for a lighter-weight image dataset. Specifically, if 𝒲_2 is obtained via, e.g., a Hungarian algorithm or a linear program, the time complexity to obtain a kernel matrix for N samples, each defined on a graph of n nodes, is roughly O(n^3 N^2). A Beckmann kernel matrix, on the other hand, runs O(n^3 + (n^2 + n)N^2) since only a single SVD calculation of L^† is required at the beginning, and then subsequent entries of the kernel matrix are obtained through matrix-vector multiplication and pairwise comparison. The Wasserstein metric 𝒲_2(α,β) was directly evaluated for each pair of examples α,β by using the CVXPY convex optimization Python package <cit.> with initialization determined by the trivial coupling βα^T. <ref> illustrates the pre-processing pipeline for an example image of a handwritten zero. 
	
	The first experiment, described in detail in <ref>, is a comparative study of ℬ_2 and 𝒲_2 in which we show that the two metrics are well-correlated for the dataset in question. The second experiment, described in <ref>, is a comparative study of ℬ_2 and 𝒲_2 within the context of unsupervised learning on the digit dataset. In this setup, we show that the two kernels perform almost identically on the dataset with respect to several metrics (with a slight edge toward ℬ_2), and thus there is some empirical evidence for drop-in usage of ℬ_2 where 𝒲_2, or perhaps 𝒲_p more generally, may present a computational bottleneck.
 
	

§ ACKNOWLEDGEMENTS


	The authors would like to acknowledge Yusu Wang for her inspirational discussions early on in this project, Stefan Steinerberger for his help placing the preliminary results within the context of the broader OT space, Zachary Lubberts for his helpful discussions on applications of graph OT to discrete geometry, and Caroline Moosmüller for her helpful feedback on linearized Wasserstein distance and Sobolev transport in the continuous setting. SR wishes to acknowledge Evangelos Nikitopoulos, as well as financial support from the Halıcıoğlu Data Science Institute through their Graduate Prize Fellowship.  AC was funded by NSF DMS 2012266 and a gift from Intel. ZW was supported by NSF CCF 2217058.
		
	plainnat
    93
    
    
    urlstyle[Abrishami et al.(2020)Abrishami, Guillen, Rule, Schutzman, Solomon, Weighill, and Wu]abrishami2020geometry
    Tara Abrishami, Nestor Guillen, Parker Rule, Zachary Schutzman, Justin Solomon, Thomas Weighill, and Si Wu.
    Geometry of graph partitions via optimal transport.
    SIAM Journal on Scientific Computing, 420 (5), 2020.
    
    [Alamgir and Luxburg(2011)]alamgir2011phase
    Morteza Alamgir and Ulrike Luxburg.
    Phase transition in the family of p-resistances.
    Advances in neural information processing systems, 24, 2011.
    
    [Alpaydin and Alimoglu(1998)]alpaydin1998penbased
    E. Alpaydin and Fevzi. Alimoglu.
    Pen-Based Recognition of Handwritten Digits.
    UCI Machine Learning Repository, 1998.
    
    [Bauer et al.(2011)Bauer, Jost, and Liu]bauer2011ollivier
    Frank Bauer, Jürgen Jost, and Shiping Liu.
    Ollivier-ricci curvature and the spectrum of the normalized graph laplace operator.
    arXiv preprint arXiv:1105.3803, 2011.
    
    [Beckmann(1952)]beckmann1952continuous
    Martin Beckmann.
    A continuous model of transportation.
    Econometrica: Journal of the Econometric Society, 1952.
    
    [Ben-Israel and Charnes(1963)]ben1963contributions
    Adi Ben-Israel and A Charnes.
    Contributions to the theory of generalized inverses.
    Journal of the Society for Industrial and Applied Mathematics, 110 (3), 1963.
    
    [Benamou and Brenier(2000)]benamou2000computational
    Jean-David Benamou and Yann Brenier.
    A computational fluid mechanics solution to the monge-kantorovich mass transfer problem.
    Numerische Mathematik, 840 (3):0 375–393, 2000.
    
    [Beveridge(2016)]beveridge2016hitting
    Andrew Beveridge.
    A hitting time formula for the discrete green's function.
    Combinatorics, Probability and Computing, 250 (3), 2016.
    
    [Beveridge and Lovász(2010)]beveridge2010exit
    Andrew Beveridge and László Lovász.
    Exit frequency matrices for finite markov chains.
    Combinatorics, Probability and Computing, 190 (4), 2010.
    
    [Bhattacharjee et al.(2023)Bhattacharjee, Cloninger, Freund, and Oslandsbotn]bhattacharjee2023effective
    Robi Bhattacharjee, Alexander Cloninger, Yoav Freund, and Andreas Oslandsbotn.
    Effective resistance in metric spaces.
    arXiv preprint arXiv:2306.15649, 2023.
    
    [Bigot(2020)]bigot2020statistical
    Jérémie Bigot.
    Statistical data analysis in the wasserstein space.
    ESAIM: Proceedings and Surveys, 68, 2020.
    
    [Black et al.(2023)Black, Wan, Nayyeri, and Wang]black2023understanding
    Mitchell Black, Zhengchao Wan, Amir Nayyeri, and Yusu Wang.
    Understanding oversquashing in gnns through the lens of effective resistance.
    In International Conference on Machine Learning, pages 2528–2547. PMLR, 2023.
    
    [Bruggner et al.(2014)Bruggner, Bodenmiller, Dill, Tibshirani, and Nolan]bruggner2014automated
    Robert V Bruggner, Bernd Bodenmiller, David L Dill, Robert J Tibshirani, and Garry P Nolan.
    Automated identification of stratifying signatures in cellular subpopulations.
    Proceedings of the National Academy of Sciences, 1110 (26), 2014.
    
    [Buitinck et al.(2013)Buitinck, Louppe, Blondel, Pedregosa, Mueller, Grisel, Niculae, Prettenhofer, Gramfort, Grobler, Layton, VanderPlas, Joly, Holt, and Varoquaux]buitinck2013sklearn
    Lars Buitinck, Gilles Louppe, Mathieu Blondel, Fabian Pedregosa, Andreas Mueller, Olivier Grisel, Vlad Niculae, Peter Prettenhofer, Alexandre Gramfort, Jaques Grobler, Robert Layton, Jake VanderPlas, Arnaud Joly, Brian Holt, and Gaël Varoquaux.
    API design for machine learning software: experiences from the scikit-learn project.
    In ECML PKDD Workshop: Languages for Data Mining and Machine Learning, 2013.
    
    [Cannon et al.(2022)Cannon, Duchin, Randall, and Rule]cannon2022spanning
    Sarah Cannon, Moon Duchin, Dana Randall, and Parker Rule.
    Spanning tree methods for sampling graph partitions.
    arXiv preprint arXiv:2210.01401, 2022.
    
    [Chen et al.(2020)Chen, Bécigneul, Ganea, Barzilay, and Jaakkola]chen2020optimal
    Benson Chen, Gary Bécigneul, Octavian-Eugen Ganea, Regina Barzilay, and Tommi Jaakkola.
    Optimal transport graph neural networks.
    arXiv preprint arXiv:2006.04804, 2020.
    
    [Chen et al.(2022)Chen, Lim, Mémoli, Wan, and Wang]chen2022weisfeiler
    Samantha Chen, Sunhyuk Lim, Facundo Mémoli, Zhengchao Wan, and Yusu Wang.
    Weisfeiler-lehman meets gromov-wasserstein.
    In International Conference on Machine Learning, pages 3371–3416. PMLR, 2022.
    
    [Chung et al.(2014)Chung, Zhao, and Kempton]chung2014ranking
    Fan Chung, Wenbo Zhao, and Mark Kempton.
    Ranking and sparsifying a connection graph.
    Internet Mathematics, 100 (1-2), 2014.
    
    [Cloninger et al.(2023a)Cloninger, Mishne, Oslandsbotn, Robertson, Wan, and Wang]cloninger2023random
    Alex Cloninger, Gal Mishne, Andreas Oslandsbotn, Sawyer Jack Robertson, Zhengchao Wan, and Yusu Wang.
    Random walks, conductance, and resistance for the connection graph laplacian.
    arXiv preprint arXiv:2308.09690, 2023a.
    
    [Cloninger et al.(2019)Cloninger, Roy, Riley, and Krumholz]cloninger2019people
    Alexander Cloninger, Brita Roy, Carley Riley, and Harlan M Krumholz.
    People mover's distance: Class level geometry using fast pairwise data adaptive transportation costs.
    Applied and Computational Harmonic Analysis, 470 (1), 2019.
    
    [Cloninger et al.(2023b)Cloninger, Hamm, Khurana, and Moosmüller]cloninger2023linearized
    Alexander Cloninger, Keaton Hamm, Varun Khurana, and Caroline Moosmüller.
    Linearized wasserstein dimensionality reduction with approximation guarantees.
    arXiv preprint arXiv:2302.07373, 2023b.
    
    [Devriendt and Lambiotte(2022)]devriendt2022discrete
    Karel Devriendt and Renaud Lambiotte.
    Discrete curvature on graphs from the effective resistance.
    Journal of Physics: Complexity, 30 (2), 2022.
    
    [Devriendt et al.(2024)Devriendt, Ottolini, and Steinerberger]devriendt2024graph
    Karel Devriendt, Andrea Ottolini, and Stefan Steinerberger.
    Graph curvature via resistance distance.
    Discrete Applied Mathematics, 348, 2024.
    
    [Diamond and Boyd(2016)]diamond2016cvxpy
    Steven Diamond and Stephen Boyd.
    CVXPY: A Python-embedded modeling language for convex optimization.
    Journal of Machine Learning Research, 2016.
    
    [Doyle and Snell(1984)]doyle1984random
    Peter G Doyle and J Laurie Snell.
    Random walks and electric networks, volume 22.
    American Mathematical Soc., 1984.
    
    [Elmoataz et al.(2008)Elmoataz, Lezoray, and Bougleux]elmoataz2008nonlocal
    Abderrahim Elmoataz, Olivier Lezoray, and Sébastien Bougleux.
    Nonlocal discrete regularization on weighted graphs: a framework for image and manifold processing.
    IEEE transactions on Image Processing, 170 (7), 2008.
    
    [Elmoataz et al.(2015)Elmoataz, Toutain, and Tenbrinck]elmoataz2015p
    Abderrahim Elmoataz, Matthieu Toutain, and Daniel Tenbrinck.
    On the p-laplacian and infinity-laplacian on graphs with applications in image and data processing.
    SIAM Journal on Imaging Sciences, 80 (4), 2015.
    
    [Essid and Solomon(2018)]essid2018quadratically
    Montacer Essid and Justin Solomon.
    Quadratically regularized optimal transport on graphs.
    SIAM Journal on Scientific Computing, 400 (4), 2018.
    
    [Evans(1997)]evans1997partial
    Lawrence C Evans.
    Partial differential equations and monge-kantorovich mass transfer.
    Current developments in mathematics, 19970 (1), 1997.
    
    [Evans and Gangbo(1999)]evans1999differential
    Lawrence C Evans and Wilfrid Gangbo.
    Differential equations methods for the Monge-Kantorovich mass transfer problem.
    American Mathematical Soc., 1999.
    
    [Fang et al.(2023)Fang, Huang, Su, and Kasai]fang2023wasserstein
    Zhongxi Fang, Jianming Huang, Xun Su, and Hiroyuki Kasai.
    Wasserstein graph distance based on l1–approximated tree edit distance between weisfeiler–lehman subtrees.
    In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, 2023.
    
    [Galichon(2018)]galichon2018optimal
    Alfred Galichon.
    Optimal transport methods in economics.
    Princeton University Press, 2018.
    
    [Gangbo and McCann(1996)]gangbo1996geometry
    Wilfrid Gangbo and Robert J McCann.
    The geometry of optimal transportation.
    Acta Mathematica, 177, 1996.
    
    [Ge and Dong(2020)]ge2020spanning
    Jun Ge and Fengming Dong.
    Spanning trees in complete bipartite graphs and resistance distance in nearly complete bipartite graphs.
    Discrete Applied Mathematics, 283, 2020.
    
    [Gerolin(2013)]arrowhead2013optimal
    Augusto Gerolin.
    Benamou-brenier’s approach for ott.
    In Optimal Transport and Applications. Summer School, Lake Arrowhead, 2013.
    
    [Gervacio(2016)]gervacio2016resistance
    Severino V Gervacio.
    Resistance distance in complete n-partite graphs.
    Discrete Applied Mathematics, 203, 2016.
    
    [Greengard et al.(2022)Greengard, Hoskins, Marshall, and Singer]greengard2022linearization
    Philip Greengard, Jeremy G Hoskins, Nicholas F Marshall, and Amit Singer.
    On a linearization of quadratic wasserstein distance.
    arXiv preprint arXiv:2201.13386, 2022.
    
    [Haker et al.(2004)Haker, Zhu, Tannenbaum, and Angenent]haker2004optimal
    Steven Haker, Lei Zhu, Allen Tannenbaum, and Sigurd Angenent.
    Optimal mass transport for registration and warping.
    International Journal of computer vision, 60, 2004.
    
    [Hubert and Arabie(1985)]hubert1985comparing
    Lawrence Hubert and Phipps Arabie.
    Comparing partitions.
    Journal of classification, 2, 1985.
    
    [Jorgensen and Pearse(2023)]jorgensen2023operator
    Palle Jorgensen and Erin PJ Pearse.
    Operator theory and analysis of infinite networks, volume 7.
    World Scientific, 2023.
    
    [Kantorovich(1942)]kantorovich1942translocation
    Leonid V Kantorovich.
    On the translocation of masses.
    In Dokl. Akad. Nauk. USSR (NS), volume 37, 1942.
    
    [Khoa and Chawla(2012)]khoa2012large
    Nguyen Lu Dang Khoa and Sanjay Chawla.
    Large scale spectral clustering using resistance distance and spielman-teng solvers.
    In Discovery Science: 15th International Conference, DS 2012, Lyon, France, October 29-31, 2012. Proceedings 15. Springer, 2012.
    
    [Khurana et al.(2023)Khurana, Kannan, Cloninger, and Moosmüller]khurana2023supervised
    Varun Khurana, Harish Kannan, Alexander Cloninger, and Caroline Moosmüller.
    Supervised learning of sheared distributions using linearized optimal transport.
    Sampling Theory, Signal Processing, and Data Analysis, 210 (1), 2023.
    
    [Kusner et al.(2015)Kusner, Sun, Kolkin, and Weinberger]kusner2015word
    Matt Kusner, Yu Sun, Nicholas Kolkin, and Kilian Weinberger.
    From word embeddings to document distances.
    In International conference on machine learning. PMLR, 2015.
    
    [Le et al.(2019)Le, Yamada, Fukumizu, and Cuturi]le2019tree
    Tam Le, Makoto Yamada, Kenji Fukumizu, and Marco Cuturi.
    Tree-sliced variants of wasserstein distances.
    Advances in neural information processing systems, 32, 2019.
    
    [Le et al.(2022)Le, Nguyen, Phung, and Nguyen]le2022sobolev
    Tam Le, Truyen Nguyen, Dinh Phung, and Viet Anh Nguyen.
    Sobolev transport: A scalable metric for probability measures with graph metrics.
    In International Conference on Artificial Intelligence and Statistics, pages 9844–9868. PMLR, 2022.
    
    [Le et al.(2023)Le, Nguyen, and Fukumizu]le2023scalable
    Tam Le, Truyen Nguyen, and Kenji Fukumizu.
    Scalable unbalanced sobolev transport for measures on a graph.
    In International Conference on Artificial Intelligence and Statistics. PMLR, 2023.
    
    [Le et al.(2024)Le, Nguyen, and Fukumizu]le2024generalized
    Tam Le, Truyen Nguyen, and Kenji Fukumizu.
    Generalized sobolev transport for probability measures on a graph.
    arXiv preprint arXiv:2402.04516, 2024.
    
    [Ling and Okada(2007)]ling2007efficient
    Haibin Ling and Kazunori Okada.
    An efficient earth mover's distance algorithm for robust histogram comparison.
    IEEE transactions on pattern analysis and machine intelligence, 290 (5), 2007.
    
    [Lovász(1993)]lovasz1993random
    László Lovász.
    Random walks on graphs.
    Combinatorics, Paul erdos is eighty, 20 (1-46), 1993.
    
    [Lovász and Winkler(1995a)]lovasz1995efficient
    László Lovász and Peter Winkler.
    Efficient stopping rules for markov chains.
    In Proceedings of the twenty-seventh annual ACM symposium on Theory of computing, 1995a.
    
    [Lovász and Winkler(1995b)]lovasz1995mixing
    László Lovász and Peter Winkler.
    Mixing of random walks and other diffusions on a graph.
    London Mathematical Society Lecture Note Series, 1995b.
    
    [Mathey-Prevot and Valette(2023)]mathey2023wasserstein
    Maxime Mathey-Prevot and Alain Valette.
    Wasserstein distance and metric trees.
    L’Enseignement Mathématique, 690 (3), 2023.
    
    [Mérigot et al.(2020)Mérigot, Delalande, and Chazal]merigot2020quantitative
    Quentin Mérigot, Alex Delalande, and Frederic Chazal.
    Quantitative stability of optimal transport maps and linearization of the 2-wasserstein space.
    In International Conference on Artificial Intelligence and Statistics. PMLR, 2020.
    
    [Miller(2007)]miller2007problem
    Stacy Miller.
    The problem of redistricting: the use of centroidal voronoi diagrams to build unbiased congressional districts.
    Senior project, Whitman College, 2007.
    
    [Monge(1781)]monge1781memoire
    Gaspard Monge.
    Mémoire sur la théorie des déblais et des remblais.
    Mem. Math. Phys. Acad. Royale Sci., 1781.
    
    [Montrucchio and Pistone(2019)]montrucchio2019kantorovich
    LUIGI Montrucchio and Giovanni Pistone.
    Kantorovich distance on a weighted graph.
    arXiv preprint arXiv:1905.07547, 1420, 2019.
    
    [Moosmüller and Cloninger(2020)]moosmuller2020linear
    Caroline Moosmüller and Alexander Cloninger.
    Linear optimal transport embedding: Provable wasserstein classification for certain rigid transformations and perturbations.
    arXiv preprint arXiv:2008.09165, 2020.
    
    [Münch and Wojciechowski(2019)]munch2019ollivier
    Florentin Münch and Radosław K Wojciechowski.
    Ollivier ricci curvature for general graph laplacians: heat equation, laplacian comparison, non-explosion and diameter bounds.
    Advances in Mathematics, 356, 2019.
    
    [Nguyen and Mamitsuka(2016)]nguyen2016new
    Canh Hao Nguyen and Hiroshi Mamitsuka.
    New resistance distances with global information on large graphs.
    In Artificial intelligence and statistics. PMLR, 2016.
    
    [Panaretos and Zemel(2019)]panaretos2019statistical
    Victor M Panaretos and Yoav Zemel.
    Statistical aspects of wasserstein distances.
    Annual review of statistics and its application, 6, 2019.
    
    [Peyré et al.(2017)Peyré, Cuturi, et al.]peyre2017computational
    Gabriel Peyré, Marco Cuturi, et al.
    Computational optimal transport.
    Center for Research in Economics and Statistics Working Papers, 10 (2017-86), 2017.
    
    [Peyre(2018)]peyre2018comparison
    Rémi Peyre.
    Comparison between w2 distance and h- 1 norm, and localization of wasserstein distance.
    ESAIM: Control, Optimisation and Calculus of Variations, 240 (4), 2018.
    
    [Ramdas et al.(2017)Ramdas, García Trillos, and Cuturi]ramdas2017wasserstein
    Aaditya Ramdas, Nicolás García Trillos, and Marco Cuturi.
    On wasserstein two-sample testing and related families of nonparametric tests.
    Entropy, 190 (2), 2017.
    
    [Robertson et al.(2023)Robertson, Kohli, Mishne, and Cloninger]robertson2023generalization
    Sawyer Robertson, Dhruv Kohli, Gal Mishne, and Alexander Cloninger.
    On a generalization of wasserstein distance and the beckmann problem to connection graphs.
    arXiv preprint arXiv:2312.10295, 2023.
    
    [Rosenberg and Hirschberg(2007)]rosenberg2007v
    Andrew Rosenberg and Julia Hirschberg.
    V-measure: A conditional entropy-based external cluster evaluation measure.
    In Proceedings of the 2007 joint conference on empirical methods in natural language processing and computational natural language learning (EMNLP-CoNLL), 2007.
    
    [Ryu et al.(2018)Ryu, Chen, Li, and Osher]ryu2018vector
    Ernest K Ryu, Yongxin Chen, Wuchen Li, and Stanley Osher.
    Vector and matrix optimal mass transport: theory, algorithm, and applications.
    SIAM Journal on Scientific Computing, 400 (5), 2018.
    
    [Saito and Herbster(2023)]saito2023multi
    Shota Saito and Mark Herbster.
    Multi-class graph clustering via approximated effective p-resistance.
    In International Conference on Machine Learning. PMLR, 2023.
    
    [Salimans et al.(2018)Salimans, Zhang, Radford, and Metaxas]salimans2018improving
    Tim Salimans, Han Zhang, Alec Radford, and Dimitris Metaxas.
    Improving gans using optimal transport.
    arXiv preprint arXiv:1803.05573, 2018.
    
    [Samal et al.(2018)Samal, Sreejith, Gu, Liu, Saucan, and Jost]samal2018comparative
    Areejit Samal, RP Sreejith, Jiao Gu, Shiping Liu, Emil Saucan, and Jürgen Jost.
    Comparative analysis of two discretizations of ricci curvature for complex networks.
    Scientific reports, 80 (1), 2018.
    
    [Santambrogio(2015)]santambrogio2015optimal
    Filippo Santambrogio.
    Optimal transport for applied mathematicians.
    Birkäuser, NY, 550 (58-63), 2015.
    
    [Schiebinger et al.(2019)Schiebinger, Shu, Tabaka, Cleary, Subramanian, Solomon, Gould, Liu, Lin, Berube, et al.]schiebinger2019optimal
    Geoffrey Schiebinger, Jian Shu, Marcin Tabaka, Brian Cleary, Vidya Subramanian, Aryeh Solomon, Joshua Gould, Siyan Liu, Stacie Lin, Peter Berube, et al.
    Optimal-transport analysis of single-cell gene expression identifies developmental trajectories in reprogramming.
    Cell, 1760 (4), 2019.
    
    [Sion(1958)]sion1958general
    Maurice Sion.
    On general minimax theorems.
    Pacific J. Math., 80 (4), 1958.
    
    [Solomon(2018)]solomon2018optimal
    Justin Solomon.
    Optimal transport on discrete domains.
    AMS Short Course on Discrete Differential Geometry, 2018.
    
    [Solomon et al.(2014)Solomon, Rustamov, Guibas, and Butscher]solomon2014earth
    Justin Solomon, Raif Rustamov, Leonidas Guibas, and Adrian Butscher.
    Earth mover's distances on discrete surfaces.
    ACM Transactions on Graphics (ToG), 330 (4), 2014.
    
    [Spielman and Srivastava(2008)]spielman2008graph
    Daniel A Spielman and Nikhil Srivastava.
    Graph sparsification by effective resistances.
    In Proceedings of the fortieth annual ACM symposium on Theory of computing, 2008.
    
    [Steinerberger(2021)]steinerberger2021wasserstein
    Stefan Steinerberger.
    A wasserstein inequality and minimal green energy on compact manifolds.
    Journal of Functional Analysis, 2810 (5), 2021.
    
    [Sugiyama and Sato(2023)]sugiyama2023kron
    Tomohiro Sugiyama and Kazuhiro Sato.
    Kron reduction and effective resistance of directed graphs.
    SIAM Journal on Matrix Analysis and Applications, 440 (1), 2023.
    
    [Takezawa et al.(2021)Takezawa, Sato, and Yamada]takezawa2021supervised
    Yuki Takezawa, Ryoma Sato, and Makoto Yamada.
    Supervised tree-wasserstein distance.
    In International Conference on Machine Learning. PMLR, 2021.
    
    [Takezawa et al.(2022)Takezawa, Sato, Kozareva, Ravi, and Yamada]takezawa2022fixed
    Yuki Takezawa, Ryoma Sato, Zornitsa Kozareva, Sujith Ravi, and Makoto Yamada.
    Fixed support tree-sliced wasserstein barycenter.
    IEICE Technical Report; IEICE Tech. Rep., 1220 (325), 2022.
    
    [TANG et al.(2017)TANG, ATHREYA, SUSSMAN, LYZINSKI, and PRIEBE]tang2017nonparametric
    MINH TANG, AVANTI ATHREYA, DANIEL L SUSSMAN, VINCE LYZINSKI, and CAREY E PRIEBE.
    A nonparametric two-sample hypothesis testing problem for random graphs.
    Bernoulli, 2017.
    
    [Tetali(1991)]tetali1991random
    Prasad Tetali.
    Random walks and the effective resistance of networks.
    Journal of Theoretical Probability, 4, 1991.
    
    [Tian et al.(2023)Tian, Lubberts, and Weber]tian2023curvature
    Yu Tian, Zachary Lubberts, and Melanie Weber.
    Curvature-based clustering on graphs.
    arXiv preprint arXiv:2307.10155, 2023.
    
    [Villani(2021)]villani2021topics
    Cédric Villani.
    Topics in optimal transportation, volume 58.
    American Mathematical Soc., 2021.
    
    [Villani et al.(2009)]villani2009optimal
    Cédric Villani et al.
    Optimal transport: old and new, volume 338.
    Springer, 2009.
    
    [Vinh et al.(2009)Vinh, Epps, and Bailey]vinh2009information
    NX Vinh, J Epps, and J Bailey.
    Information theoretic measures for clusterings comparison: Variants.
    Properties, Normalization and Correction for Chance, 18, 2009.
    
    [Von Luxburg(2007)]von2007tutorial
    Ulrike Von Luxburg.
    A tutorial on spectral clustering.
    Statistics and computing, 17, 2007.
    
    [Wang et al.(2013)Wang, Slepčev, Basu, Ozolek, and Rohde]wang2013linear
    Wei Wang, Dejan Slepčev, Saurav Basu, John A Ozolek, and Gustavo K Rohde.
    A linear optimal transportation framework for quantifying and visualizing variations in sets of images.
    International journal of computer vision, 101, 2013.
    
    [Weaver(2018)]weaver2018lipschitz
    Nik Weaver.
    Lipschitz algebras.
    World Scientific, 2018.
    
    [Yamada et al.(2022)Yamada, Takezawa, Sato, Bao, Kozareva, and Ravi]yamada2022approximating
    Makoto Yamada, Yuki Takezawa, Ryoma Sato, Han Bao, Zornitsa Kozareva, and Sujith Ravi.
    Approximating 1-wasserstein distance with trees.
    Transactions on Machine Learning Research, 2022.
    
    [Ye and Yan(2019)]ye2019resistance
    Luzh Ye and Weigen Yan.
    Resistance between two vertices of almost complete bipartite graphs.
    Discrete Applied Mathematics, 257, 2019.
    
    [Ye et al.(2021)Ye, Ma, Ni, Liu, Gao, and Chen]ye2021riccignn
    Ze Ye, Tengfei Ma, Chien-Chun Ni, Kin Sum Liu, Jie Gao, and Chao Chen.
    Ricci-{gnn}: Defending against structural attacks through a geometric approach, 2021.
    URL <https://openreview.net/forum?id=_qoQkWNEhS>.
    
    [Zhang et al.(2010)Zhang, Jin, and Zhou]zhang2010understanding
    Yin Zhang, Rong Jin, and Zhi-Hua Zhou.
    Understanding bag-of-words model: a statistical framework.
    International journal of machine learning and cybernetics, 1, 2010.
    
    

 
 





§ PROOFS


This appendix contains deferred proofs from the remainder of the paper.



 §.§ Proofs from section 2



	We can derive the Lagrangian dual in a fairly straightforward manner. We rewrite the constraint BJ = α-β using the expresion
	
    sup_ϕ∈ℓ(V)ϕ^T(α-β) -  ϕ^T (BJ)    =
    			0     if BJ = α-β
    
    			+∞    otherwise.
 
	Thus the Beckmann problem may be rewritten as the following
	
    inf_BJ = α-βJ_w, p   = inf_J {J_w, p + sup_ϕϕ^T(α-β) - ϕ^T(BJ) }
       =sup_ϕ{ϕ^T(α-β) + inf_J J_w, p -ϕ^T(BJ)},

	where we exchange the sup and inf using, e.g., Sion's minimax theorem <cit.> or Slater's condition. Now we have that, using the Legendre transform of ·_p,
	
    inf_J  J(e)_w,p - ϕ^T(BJ)    = inf_J  J_w, p - J^T (B^Tϕ)
       = 
    			0     if B^Tϕ_w^1-q,q≤ 1
    
    			-∞    otherwise

	since the dual norm to ·_w, p is ·_w^1-q, q. Therefore, the Lagrangian dual becomes
	
    sup{ϕ^T(α-β) : ϕ∈ℓ(V), B^Tϕ_w^1-q, q≤ 1 }

	The special case of p=1 and q=∞ follows from the same proof, with the only change being that the dual norm to ·_w, 1 is ·_w^-1, ∞, which can be shown directly or through a limiting argument.



	Let i, j∈ V be fixed distinct vertices. We begin by noting simply that by inspection,
	
    B(π_ij I_P_ij) = π_ij(δ_i - δ_j)

	and thus that
	
    B(∑_i, j∈ Vπ_ij I_P_ij)    = ∑_i,j ∈ Vπ_ij (δ_i - δ_j)
       = ∑_i∈ Vδ_i (π1)_i - ∑_j∈ Vδ_j (1^T π)_j
       = α - β.

	Since J^π≥ 0 by construction, J^π is a feasible edge flow for ℬ_p(α,β). The estimate on J^π_w,p follows immediately from the triangle inequality.



	Suppose π∈Π(α,β) is an optimal coupling for 𝒲_d_p, p(α,β). For each i, j, let P_ij∈𝒫(i, j) be a choice of path which is minimal in the sense of d_p, and let J^π = J^π, (P_ij) be the feasible flow as in <ref>. Then by Hölder's inequality,
	
    ℬ_p(α,β)    ≤J^π_w,p
       ≤∑_i, j∈ Vπ_ijI_P_ij_w, p
       ≤∑_i, j∈ Vπ_ij d_p(i, j) 
       ≤1_V× V_q(∑_i, j∈ Vπ_ij^p d_p(i, j)^p)^1/p
       ≤1_V× V_q(∑_i, j∈ Vπ_ij d_p(i, j)^p)^1/p = n^2/q𝒲_d_p, p(α,β)

	since π_ij≤ 1 and p >1.



    The proof of this result is straightforward and comes in two parts. First, we observe that 𝒲_p(α,β) = 𝒲_p(α, β) where
        
    𝒲_p(α, β)^p = inf{∫_ℝ×ℝ |x-y|^p π(x, y) : π∈Π(α, β) }.
 
    In this case, Π(α, β) is the set of all Borel probability measures on ℝ×ℝ with respective marginals α,β. Equality here holds because any such π∈Π(α, β), since it has discrete marginals, must itself be discrete and supported on the Cartesian product {1, …, n}×{1, …, n}; i.e., it can be written as a matrix. Thus, there is a one-to-one correspondence between the feasible couplings on V× V and ℝ×ℝ, and the cost of any such pair of corresponding couplings is identical since shortest path distance on V(P_n) is |i-j|.
    
    The second part of the proof is classical, and it suffices to observe that the p-Wasserstein distance between probability measures on ℝ is simply the p-norm distance between their inverse cdfs. For a proof we recommend, e.g., <cit.>.



    We note that the signed incidence matrix B, which is of shape n× n-1, must have no kernel- for the rank of B agrees with the rank of BB^T = L, which is n-1 for any connected graph on n nodes. Therefore, for any 1≤ p<∞, ℬ_p(α,β) is the p-norm of the unique flow J such that BJ = α-β and which does not depend on p. We therefore need only show that B(K_α - K_β) =α-β. To wit, if i = 1, note that
        
    BK_α(1) = K_α(1, 2) = α(1),

    and if 1 < i <n, 
        
    BK_α(i) = K_α(i, i+1) - K_α(i-1, i) = α(i),

    and lastly if i = n,
        
    BK_α(n) = -K_α(n-1, n) = α(n) - 1.

    Therefore, with a simple cancellation, it holds B(K_α - K_β) =α-β for all i∈ V. The first equality in the proposition therefore follows, and the second follows immediately upon inspection.





    Once again the proof of this result is a matter of rigid feasibility: there can only really be one feasible flow J∈ℓ(E') satisfying BJ = α-β owing to, e.g., rank considerations for Laplacians on trees. Thus it is enough only to establish that B(K_α-K_β) = α-β. Let i∈ V be fixed, and suppose i has N incoming oriented edges, and M outgoing oriented edges. At the root of each incoming oriented edge is a subtree I_ℓ, and at the head of each outgoing oriented edge is a subtree O_ℓ. This setup is illustrated in <ref>. For each 1≤ℓ≤ N, let α(I_ℓ) be the total mass of α on that component, and similarly for α(O_ℓ). We split up the argument into a few cases depending on M, N. Assume for a moment that M≥ 2 and N≥ 1; then by inspection we have
        
    (BK_α)(i)    = -(α(I_1) + α(I_2) + … + α(I_N)) 
       + (α(i) + ∑_ℓ=1^N α(I_ℓ) + ∑_ℓ≠ 1α(O_ℓ) ) + ⋯ + (α(i) + ∑_ℓ=1^N α(I_ℓ) + ∑_ℓ≠ Mα(O_ℓ) )
       = Mα(i) + (M-1)∑_ℓ=1^N α(I_ℓ) + (M-1) ∑_ℓ=1^M α(O_ℓ) 
       = α(i) + M - 1.

    If M=1 and N≥ 1, then (BK_α)(i) = α(i). If M = 0 and N≥ 1, then (BK_α)(i) = α(i) - 1. If N=0, then (BK_α)(i) = α(i) + M - 1, as before, since the only difference is the lack of contribution of the I_ℓ's. Therefore, BK_α = α + c where c∈ℝ^n is a vector which does not depend on α, and thus it holds that B(K_α-K_β) = α-β. The claim follows.




 §.§ Proofs from section 3



	We start by considering the pointwise value of the vector L^† (α-β). To this end, a useful observation is that for any coupling π∈Π(α,β), we have
	
    α-β = ∑_i, j∈ Vπ_ij (δ_i - δ_j).

	With this in hand alongside <ref>, we can calculate
	
    L^† (α-β)_i    = ∑_j∈ V( (L^†)_ij∑_k, ℓ∈ Vπ_kℓ(δ_k(j) - δ_ℓ(j)) )
       =1/(G)∑_j, k, ℓ∈ V( H(ρ, j) - H(i, j)) π_kℓ(δ_k(j) - δ_ℓ(j))
       = 1/(G)∑_k, ℓ∈ Vπ_kℓ( H(ρ, k) - H(i, k) - H(ρ, ℓ) + H(i, ℓ)  )
       = 1/(G)∑_kα_k (H(ρ, k) - H(i, k)) - ∑_ℓβ_ℓ ( H(ρ, ℓ) - H(i,ℓ) ) 
       = 1/(G)∑_k (α_k - β_k) (H(ρ, k) - H(i, k)) 
       = c - 1/(G)∑_k (α_k - β_k) H(i, k)

	where by introducing the placeholder c for some c∈ℝ, we are separating the first piece of the summand which does not depend on i and which will be eliminated in the subsequent calculation. Then, multiplying by (α-β)^T, and using the fact that (α-β)^T 1_n = 0, we get
    
	
    (α-β)^T L^† (α-β)    = -1/(G)∑_i, k∈ V(α_i - β_i)   (α_k - β_k) H(i, k) 
       = -1/(G)∑_k∈ V (α_k - β_k) ∑_i∈ V( α_i H(i, k) - β_i H(i, k) )
       =  -1/(G)∑_k∈ V (α_k - β_k)(H(α, k) - H(β, k))

	using <ref>(i) in the final line. The claim follows.



	From the proof of <ref> and Hölder's inequality, we have
	
    r_αβ   =-1/(G)∑_i, k∈ V(α_i - β_i)   (α_k - β_k) H(i, k) 
       =1/(G)∑_k   (β_k - α_k) ∑_i  (α_i - β_i)  H(i, k)
       ≤1/ (G)α-β_1 ∑_i  (α_i - β_i)  H(i, k)_∞
       ≤2/(G)max{H(α,β), H(β,α)}.

    
	where in the final line we make use of <ref>(ii). The second inequality in the Corollary follows from the observation in the proof of <ref> that if π∈Π(α,β) is any coupling, α-β = ∑_i, j∈ Vπ_ij (δ_i - δ_j). Then since x↦L^-1/2 x_2^2 is convex in x, 
	
    (α-β)^T L^† (α-β)    = L^-1/2 (α-β)_2^2
       = ∑_i, jπ_ij L^-1/2(δ_i-δ_j)_2^2
       ≤∑_i, jπ_ij r_ij.

	Thus in the special case where π = αβ^T,
	
    r_αβ   ≤1/(G)∑_i, j∈ Vα_iβ_j  (H(i, j) + H(j, i))
        = 1/(G) (H_n(α,β) + H_n(β, α)).



 
	Let μ_t ∈ C^1([0, 1]) such that μ_0 =α and μ_1 = β. Then we estimate, using <ref> and the fact that L^† = L^† L L^†,
	
    ∫_0^1 μ_t _Ḣ^-1(V)^2  t    = ∫_0^1 (μ_t)^T L^†μ_t  t
       =  ∫_0^1 (L^†μ_t)^T L (L^†μ_t)  t
       = ∫_0^1 ∑_(i, j)∈ E' w_ij|L^†μ_t(i) - L^†μ_t(j)|^2  t
       = ∑_(i, j)∈ E' w_ij∫_0^1 |L^†μ_t(i) - L^†μ_t(j)|^2 t
       ≥∑_(i, j)∈ E' w_ij| ∫_0^1 (L^†μ_t(i) - L^†μ_t(j))  t |^2
 
	where the final inequality follows from Jensen's inequality. Then, focusing on the inner term and letting e=(i, j)∈ E' we evaluate
	
    ∫_0^1 (L^†μ_t(i) - L^†μ_t(j))  t    = B^T(∫_0^1 L^†μ_t  t)(e)
       =B^T L^†(∫_0^1 μ_t  t)(e)
       = B^T L^†(β-α)(e).

	Hence,
	
    ∑_(i, j)∈ E' w_ij| ∫_0^1 (L^†μ_t(i) - L^†μ_t(j))  t |^2     = ∑_(i, j)∈ E' w_ij|B^T L^†(β-α)(i, j)|^2
       = (L^† (β-α))^T L  (L^† (β-α))
       =(α-β)^T L^† (α-β),

	and therefore,
	
    inf{∫_0^1 μ_t _Ḣ^-1(V)^2  t: μ_t∈ C^1([0, 1]), μ_0 = α, μ_1 = β}≥ (α-β)^T L^† (α-β) = ℬ_2(α,β)^2.

	For the reverse inequality, it is sufficient to consider the special curve μ_t = (1-t)α - tβ, which satisfies ∫_0^1 μ_t _Ḣ^-1(V)^2  t = (α-β)^T L^† (α-β) = ℬ_2(α,β)^2.



	We observe first that the collection of vectors { T}_T∈𝒯⊂ℓ(V) is convex; for the first two conditions (i) and (ii) are equivalent to the intersection of a finite number of linear constraints; and (iii) is stable under convex combinations. Since L^-1/2 is a linear map, the image L^-1/2𝒜_i will be convex for i=1, 2. Separately, we have for S, T∈𝒯:
	
    Tα_1 - Sα_2_2    =  (α_1 - α_2) + ( T -  S)_2 
       ≥| α_1 - α_2_2 -  T -  S_2 |
       ≥α_1 - α_2_2 - 2δ > 0.

	Thus, 
	
    min_S, T∈𝒯Tα_1 - Sα_2_2 >0

	and hence 𝒜_1∩𝒜_2  = ∅. Note moreover that L^-1/2 (and more generally, any power of L) will act nonsingularly on mean zero vectors since its kernel is exactly the constant functions. Therefore, since Tα_1 - Sα_2 is mean zero for any T, S, and since 𝒜_1 and 𝒜_2 are disjoint, the sets L^-1/2𝒜_1 and L^-1/2𝒜_2 are disjoint as well; and thus they are linearly separable.




















































