Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly

Sign up
Appearance settings

Commit 55eecd1

Browse files
Introducing Sinkhorn's method, also started to improve LP explanation for Optimal Transport
1 parent aeed5d9 commit 55eecd1

File tree

1 file changed

+56
-3
lines changed

1 file changed

+56
-3
lines changed

‎OptimalTransportWasserteinDistance.ipynb‎

Lines changed: 56 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -225,7 +225,18 @@
225225
" \\mathrm{so \\ that} \\ & \\mathbf{A}^T \\mathbf{y} & \\leq & \\ \\mathbf{c} \\\\ \\\\\n",
226226
" \\end{array}\n",
227227
"\\end{array}\n",
228-
"$$"
228+
"$$\n",
229+
"The fact that $\\tilde{z}$ is a lower bound for $z$ (weak duality) immediately comes to play here:\n",
230+
"$$\n",
231+
"\\begin{align*}\n",
232+
" c & \\geq A^T y \\\\\n",
233+
" c^T & \\geq y^T A \\\\\n",
234+
" c^T x & \\geq y^T A x \\\\\n",
235+
" c^T x & \\geq y^T b \\\\\n",
236+
" z & \\geq \\tilde{z}\n",
237+
"\\end{align*}\n",
238+
"$$\n",
239+
"Strong duality is a bit more complex to prove, you'll find elements everywhere, for instance on Vincent Therrmann [blog](https://vincentherrmann.github.io/blog/wasserstein/) "
229240
]
230241
},
231242
{
@@ -286,9 +297,51 @@
286297
"cell_type": "markdown",
287298
"metadata": {},
288299
"source": [
289-
"## Probabilistic formulation, doubly stochastic matrices and Sinkhorn theorem\n",
300+
"## Probabilistic formulation, doubly stochastic matrices and Sinkhorn's theorem\n",
301+
"\n",
302+
"This part is largly inspired by the following wikipedia article on [sinkhorn theorem](https://en.wikipedia.org/wiki/Sinkhorn%27s_theorem) and [doubly stochastic matrices](https://en.wikipedia.org/wiki/Doubly_stochastic_matrix) as well as this [blog post](https://zipjiang.github.io/2020/11/23/sinkhorn's-theorem-,-sinkhorn-algorithm-and-applications.html)\n",
303+
"\n",
304+
"### (Doubly) Stochastic matrices\n",
305+
"A doubly stochastic matrix (also called bistochastic matrix) is a square matrix $X = (x_{ij})$ of nonnegative real numbers, each of whose rows and columns sums to 1, i.e: $ \\sum_{i}x_{ij}=\\sum_{j}x_{ij}=1 $\n",
306+
"Thus, a doubly stochastic matrix is both left stochastic and right stochastic.\n",
307+
"Any matrix that is both left and right stochastic must be square: if every row sums to 1 then the sum of all entries in the matrix must be equal to the number of rows, and since the same holds for columns, the number of rows and columns must be equal.\n",
308+
"\n",
309+
"### Sinkhorn's theorem\n",
310+
"If $A$ is an $n \\times n$ matrix with strictly positive elements, then there exist diagonal matrices $D_1$ and $D_2$ with strictly positive diagonal elements such that $D_1AD_2$ is doubly stochastic. The matrices $D_1$ and $D_2$ are unique modulo multiplying the first matrix by a positive number and dividing the second one by the same number.\n",
311+
"\n",
312+
"### Sinkhorn's algorithm\n",
313+
"What is great about Sinkhorn’s theorem is that there is a efficient algorithm to calculate an approximation of the doubly stochastic matrix given a matrix $A,ドル that features linear convergence. And the algorithm is very simple: one just iteratively normalize the matrix along rows and columns."
314+
]
315+
},
316+
{
317+
"cell_type": "code",
318+
"execution_count": null,
319+
"metadata": {},
320+
"outputs": [],
321+
"source": [
322+
"# Sinkhorn method\n",
323+
"\n",
324+
"def sinkhorn(A, L):\n",
325+
" # Code for calculating the doubly stochastic matrix\n",
326+
" # using Sinkhorn-Knopp algorithm.\n",
327+
" # ----------\n",
328+
" # Input: positive matrix A[N x N], max iteration L\n",
329+
"\n",
330+
" for i in range(L):\n",
331+
" A = A / np.matmul(A, np.ones(N, 1))\n",
332+
" A = A / np.matmul(np.ones(1, N), A)\n",
333+
"\n",
334+
" row_residual = np.dot(A, np.ones(N))-1\n",
335+
" col_residual = np.dot(A.T, np.ones(N))-1\n",
336+
"\n",
337+
" # Test for convergence and early stop.\n",
338+
" if np.allclose(col_residual, 0) and np.allclose(row_residual, 0)::\n",
339+
" break\n",
340+
"\n",
341+
" return A\n",
290342
"\n",
291-
"This part is largly inspired by the following wikipedia article on [sinkhorn theorem](https://en.wikipedia.org/wiki/Sinkhorn%27s_theorem) and [doubly stochastic matrices](https://en.wikipedia.org/wiki/Doubly_stochastic_matrix)\n"
343+
"N = 10\n",
344+
"A = np.random.uniform(0,100,(N,N))"
292345
]
293346
},
294347
{

0 commit comments

Comments
(0)

AltStyle によって変換されたページ (->オリジナル) /