add all pages in documentation

7 files changed, 1147 insertions, 0 deletions

diff --git a/docs/source/auto_examples/plot_otda_jcpot.ipynb b/docs/source/auto_examples/plot_otda_jcpot.ipynb
new file mode 100644
index 0000000..a81d47a
--- /dev/null
+++ b/docs/source/auto_examples/plot_otda_jcpot.ipynb
@@ -0,0 +1,173 @@
+{
+  "cells": [
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "%matplotlib inline"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "\n# OT for multi-source target shift\n\n\nThis example introduces a target shift problem with two 2D source and 1 target domain.\n\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "# Authors: Remi Flamary <remi.flamary@unice.fr>\n#          Ievgen Redko <ievgen.redko@univ-st-etienne.fr>\n#\n# License: MIT License\n\nimport pylab as pl\nimport numpy as np\nimport ot\nfrom ot.datasets import make_data_classif"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Generate data\n-------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "n = 50\nsigma = 0.3\nnp.random.seed(1985)\n\np1 = .2\ndec1 = [0, 2]\n\np2 = .9\ndec2 = [0, -2]\n\npt = .4\ndect = [4, 0]\n\nxs1, ys1 = make_data_classif('2gauss_prop', n, nz=sigma, p=p1, bias=dec1)\nxs2, ys2 = make_data_classif('2gauss_prop', n + 1, nz=sigma, p=p2, bias=dec2)\nxt, yt = make_data_classif('2gauss_prop', n, nz=sigma, p=pt, bias=dect)\n\nall_Xr = [xs1, xs2]\nall_Yr = [ys1, ys2]"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "da = 1.5\n\n\ndef plot_ax(dec, name):\n    pl.plot([dec[0], dec[0]], [dec[1] - da, dec[1] + da], 'k', alpha=0.5)\n    pl.plot([dec[0] - da, dec[0] + da], [dec[1], dec[1]], 'k', alpha=0.5)\n    pl.text(dec[0] - .5, dec[1] + 2, name)"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Fig 1 : plots source and target samples\n---------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "pl.figure(1)\npl.clf()\nplot_ax(dec1, 'Source 1')\nplot_ax(dec2, 'Source 2')\nplot_ax(dect, 'Target')\npl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9,\n           label='Source 1 ({:1.2f}, {:1.2f})'.format(1 - p1, p1))\npl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9,\n           label='Source 2 ({:1.2f}, {:1.2f})'.format(1 - p2, p2))\npl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9,\n           label='Target ({:1.2f}, {:1.2f})'.format(1 - pt, pt))\npl.title('Data')\n\npl.legend()\npl.axis('equal')\npl.axis('off')"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Instantiate Sinkhorn transport algorithm and fit them for all source domains\n----------------------------------------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "ot_sinkhorn = ot.da.SinkhornTransport(reg_e=1e-1, metric='sqeuclidean')\n\n\ndef print_G(G, xs, ys, xt):\n    for i in range(G.shape[0]):\n        for j in range(G.shape[1]):\n            if G[i, j] > 5e-4:\n                if ys[i]:\n                    c = 'b'\n                else:\n                    c = 'r'\n                pl.plot([xs[i, 0], xt[j, 0]], [xs[i, 1], xt[j, 1]], c, alpha=.2)"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Fig 2 : plot optimal couplings and transported samples\n------------------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "pl.figure(2)\npl.clf()\nplot_ax(dec1, 'Source 1')\nplot_ax(dec2, 'Source 2')\nplot_ax(dect, 'Target')\nprint_G(ot_sinkhorn.fit(Xs=xs1, Xt=xt).coupling_, xs1, ys1, xt)\nprint_G(ot_sinkhorn.fit(Xs=xs2, Xt=xt).coupling_, xs2, ys2, xt)\npl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)\npl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)\npl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)\n\npl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')\npl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')\n\npl.title('Independent OT')\n\npl.legend()\npl.axis('equal')\npl.axis('off')"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Instantiate JCPOT adaptation algorithm and fit it\n----------------------------------------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "otda = ot.da.JCPOTTransport(reg_e=1, max_iter=1000, metric='sqeuclidean', tol=1e-9, verbose=True, log=True)\notda.fit(all_Xr, all_Yr, xt)\n\nws1 = otda.proportions_.dot(otda.log_['D2'][0])\nws2 = otda.proportions_.dot(otda.log_['D2'][1])\n\npl.figure(3)\npl.clf()\nplot_ax(dec1, 'Source 1')\nplot_ax(dec2, 'Source 2')\nplot_ax(dect, 'Target')\nprint_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)\nprint_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)\npl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)\npl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)\npl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)\n\npl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')\npl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')\n\npl.title('OT with prop estimation ({:1.3f},{:1.3f})'.format(otda.proportions_[0], otda.proportions_[1]))\n\npl.legend()\npl.axis('equal')\npl.axis('off')"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Run oracle transport algorithm with known proportions\n----------------------------------------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "h_res = np.array([1 - pt, pt])\n\nws1 = h_res.dot(otda.log_['D2'][0])\nws2 = h_res.dot(otda.log_['D2'][1])\n\npl.figure(4)\npl.clf()\nplot_ax(dec1, 'Source 1')\nplot_ax(dec2, 'Source 2')\nplot_ax(dect, 'Target')\nprint_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)\nprint_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)\npl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)\npl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)\npl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)\n\npl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')\npl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')\n\npl.title('OT with known proportion ({:1.1f},{:1.1f})'.format(h_res[0], h_res[1]))\n\npl.legend()\npl.axis('equal')\npl.axis('off')\npl.show()"
+      ]
+    }
+  ],
+  "metadata": {
+    "kernelspec": {
+      "display_name": "Python 3",
+      "language": "python",
+      "name": "python3"
+    },
+    "language_info": {
+      "codemirror_mode": {
+        "name": "ipython",
+        "version": 3
+      },
+      "file_extension": ".py",
+      "mimetype": "text/x-python",
+      "name": "python",
+      "nbconvert_exporter": "python",
+      "pygments_lexer": "ipython3",
+      "version": "3.6.9"
+    }
+  },
+  "nbformat": 4,
+  "nbformat_minor": 0
+}
\ No newline at end of file
diff --git a/docs/source/auto_examples/plot_otda_jcpot.py b/docs/source/auto_examples/plot_otda_jcpot.py
new file mode 100644
index 0000000..c495690
--- /dev/null
+++ b/docs/source/auto_examples/plot_otda_jcpot.py
@@ -0,0 +1,171 @@
+# -*- coding: utf-8 -*-
+"""
+========================
+OT for multi-source target shift
+========================
+
+This example introduces a target shift problem with two 2D source and 1 target domain.
+
+"""
+
+# Authors: Remi Flamary <remi.flamary@unice.fr>
+#          Ievgen Redko <ievgen.redko@univ-st-etienne.fr>
+#
+# License: MIT License
+
+import pylab as pl
+import numpy as np
+import ot
+from ot.datasets import make_data_classif
+
+##############################################################################
+# Generate data
+# -------------
+n = 50
+sigma = 0.3
+np.random.seed(1985)
+
+p1 = .2
+dec1 = [0, 2]
+
+p2 = .9
+dec2 = [0, -2]
+
+pt = .4
+dect = [4, 0]
+
+xs1, ys1 = make_data_classif('2gauss_prop', n, nz=sigma, p=p1, bias=dec1)
+xs2, ys2 = make_data_classif('2gauss_prop', n + 1, nz=sigma, p=p2, bias=dec2)
+xt, yt = make_data_classif('2gauss_prop', n, nz=sigma, p=pt, bias=dect)
+
+all_Xr = [xs1, xs2]
+all_Yr = [ys1, ys2]
+# %%
+
+da = 1.5
+
+
+def plot_ax(dec, name):
+    pl.plot([dec[0], dec[0]], [dec[1] - da, dec[1] + da], 'k', alpha=0.5)
+    pl.plot([dec[0] - da, dec[0] + da], [dec[1], dec[1]], 'k', alpha=0.5)
+    pl.text(dec[0] - .5, dec[1] + 2, name)
+
+
+##############################################################################
+# Fig 1 : plots source and target samples
+# ---------------------------------------
+
+pl.figure(1)
+pl.clf()
+plot_ax(dec1, 'Source 1')
+plot_ax(dec2, 'Source 2')
+plot_ax(dect, 'Target')
+pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9,
+           label='Source 1 ({:1.2f}, {:1.2f})'.format(1 - p1, p1))
+pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9,
+           label='Source 2 ({:1.2f}, {:1.2f})'.format(1 - p2, p2))
+pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9,
+           label='Target ({:1.2f}, {:1.2f})'.format(1 - pt, pt))
+pl.title('Data')
+
+pl.legend()
+pl.axis('equal')
+pl.axis('off')
+
+##############################################################################
+# Instantiate Sinkhorn transport algorithm and fit them for all source domains
+# ----------------------------------------------------------------------------
+ot_sinkhorn = ot.da.SinkhornTransport(reg_e=1e-1, metric='sqeuclidean')
+
+
+def print_G(G, xs, ys, xt):
+    for i in range(G.shape[0]):
+        for j in range(G.shape[1]):
+            if G[i, j] > 5e-4:
+                if ys[i]:
+                    c = 'b'
+                else:
+                    c = 'r'
+                pl.plot([xs[i, 0], xt[j, 0]], [xs[i, 1], xt[j, 1]], c, alpha=.2)
+
+
+##############################################################################
+# Fig 2 : plot optimal couplings and transported samples
+# ------------------------------------------------------
+pl.figure(2)
+pl.clf()
+plot_ax(dec1, 'Source 1')
+plot_ax(dec2, 'Source 2')
+plot_ax(dect, 'Target')
+print_G(ot_sinkhorn.fit(Xs=xs1, Xt=xt).coupling_, xs1, ys1, xt)
+print_G(ot_sinkhorn.fit(Xs=xs2, Xt=xt).coupling_, xs2, ys2, xt)
+pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+pl.title('Independent OT')
+
+pl.legend()
+pl.axis('equal')
+pl.axis('off')
+
+##############################################################################
+# Instantiate JCPOT adaptation algorithm and fit it
+# ----------------------------------------------------------------------------
+otda = ot.da.JCPOTTransport(reg_e=1, max_iter=1000, metric='sqeuclidean', tol=1e-9, verbose=True, log=True)
+otda.fit(all_Xr, all_Yr, xt)
+
+ws1 = otda.proportions_.dot(otda.log_['D2'][0])
+ws2 = otda.proportions_.dot(otda.log_['D2'][1])
+
+pl.figure(3)
+pl.clf()
+plot_ax(dec1, 'Source 1')
+plot_ax(dec2, 'Source 2')
+plot_ax(dect, 'Target')
+print_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)
+print_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)
+pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+pl.title('OT with prop estimation ({:1.3f},{:1.3f})'.format(otda.proportions_[0], otda.proportions_[1]))
+
+pl.legend()
+pl.axis('equal')
+pl.axis('off')
+
+##############################################################################
+# Run oracle transport algorithm with known proportions
+# ----------------------------------------------------------------------------
+h_res = np.array([1 - pt, pt])
+
+ws1 = h_res.dot(otda.log_['D2'][0])
+ws2 = h_res.dot(otda.log_['D2'][1])
+
+pl.figure(4)
+pl.clf()
+plot_ax(dec1, 'Source 1')
+plot_ax(dec2, 'Source 2')
+plot_ax(dect, 'Target')
+print_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)
+print_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)
+pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+pl.title('OT with known proportion ({:1.1f},{:1.1f})'.format(h_res[0], h_res[1]))
+
+pl.legend()
+pl.axis('equal')
+pl.axis('off')
+pl.show()
diff --git a/docs/source/auto_examples/plot_otda_jcpot.rst b/docs/source/auto_examples/plot_otda_jcpot.rst
new file mode 100644
index 0000000..3433190
--- /dev/null
+++ b/docs/source/auto_examples/plot_otda_jcpot.rst
@@ -0,0 +1,336 @@
+.. only:: html
+
+    .. note::
+        :class: sphx-glr-download-link-note
+
+        Click :ref:`here <sphx_glr_download_auto_examples_plot_otda_jcpot.py>`     to download the full example code
+    .. rst-class:: sphx-glr-example-title
+
+    .. _sphx_glr_auto_examples_plot_otda_jcpot.py:
+
+
+========================
+OT for multi-source target shift
+========================
+
+This example introduces a target shift problem with two 2D source and 1 target domain.
+
+
+
+.. code-block:: default
+
+
+    # Authors: Remi Flamary <remi.flamary@unice.fr>
+    #          Ievgen Redko <ievgen.redko@univ-st-etienne.fr>
+    #
+    # License: MIT License
+
+    import pylab as pl
+    import numpy as np
+    import ot
+    from ot.datasets import make_data_classif
+
+
+
+
+
+
+
+
+Generate data
+-------------
+
+
+.. code-block:: default
+
+    n = 50
+    sigma = 0.3
+    np.random.seed(1985)
+
+    p1 = .2
+    dec1 = [0, 2]
+
+    p2 = .9
+    dec2 = [0, -2]
+
+    pt = .4
+    dect = [4, 0]
+
+    xs1, ys1 = make_data_classif('2gauss_prop', n, nz=sigma, p=p1, bias=dec1)
+    xs2, ys2 = make_data_classif('2gauss_prop', n + 1, nz=sigma, p=p2, bias=dec2)
+    xt, yt = make_data_classif('2gauss_prop', n, nz=sigma, p=pt, bias=dect)
+
+    all_Xr = [xs1, xs2]
+    all_Yr = [ys1, ys2]
+
+
+
+
+
+
+
+
+.. code-block:: default
+
+
+    da = 1.5
+
+
+    def plot_ax(dec, name):
+        pl.plot([dec[0], dec[0]], [dec[1] - da, dec[1] + da], 'k', alpha=0.5)
+        pl.plot([dec[0] - da, dec[0] + da], [dec[1], dec[1]], 'k', alpha=0.5)
+        pl.text(dec[0] - .5, dec[1] + 2, name)
+
+
+
+
+
+
+
+
+
+Fig 1 : plots source and target samples
+---------------------------------------
+
+
+.. code-block:: default
+
+
+    pl.figure(1)
+    pl.clf()
+    plot_ax(dec1, 'Source 1')
+    plot_ax(dec2, 'Source 2')
+    plot_ax(dect, 'Target')
+    pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9,
+               label='Source 1 ({:1.2f}, {:1.2f})'.format(1 - p1, p1))
+    pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9,
+               label='Source 2 ({:1.2f}, {:1.2f})'.format(1 - p2, p2))
+    pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9,
+               label='Target ({:1.2f}, {:1.2f})'.format(1 - pt, pt))
+    pl.title('Data')
+
+    pl.legend()
+    pl.axis('equal')
+    pl.axis('off')
+
+
+
+
+.. image:: /auto_examples/images/sphx_glr_plot_otda_jcpot_001.png
+    :class: sphx-glr-single-img
+
+
+.. rst-class:: sphx-glr-script-out
+
+ Out:
+
+ .. code-block:: none
+
+
+    (-1.85, 5.85, -4.1171725099266725, 4.197384527473105)
+
+
+
+Instantiate Sinkhorn transport algorithm and fit them for all source domains
+----------------------------------------------------------------------------
+
+
+.. code-block:: default
+
+    ot_sinkhorn = ot.da.SinkhornTransport(reg_e=1e-1, metric='sqeuclidean')
+
+
+    def print_G(G, xs, ys, xt):
+        for i in range(G.shape[0]):
+            for j in range(G.shape[1]):
+                if G[i, j] > 5e-4:
+                    if ys[i]:
+                        c = 'b'
+                    else:
+                        c = 'r'
+                    pl.plot([xs[i, 0], xt[j, 0]], [xs[i, 1], xt[j, 1]], c, alpha=.2)
+
+
+
+
+
+
+
+
+
+Fig 2 : plot optimal couplings and transported samples
+------------------------------------------------------
+
+
+.. code-block:: default
+
+    pl.figure(2)
+    pl.clf()
+    plot_ax(dec1, 'Source 1')
+    plot_ax(dec2, 'Source 2')
+    plot_ax(dect, 'Target')
+    print_G(ot_sinkhorn.fit(Xs=xs1, Xt=xt).coupling_, xs1, ys1, xt)
+    print_G(ot_sinkhorn.fit(Xs=xs2, Xt=xt).coupling_, xs2, ys2, xt)
+    pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+    pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+    pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+    pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+    pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+    pl.title('Independent OT')
+
+    pl.legend()
+    pl.axis('equal')
+    pl.axis('off')
+
+
+
+
+.. image:: /auto_examples/images/sphx_glr_plot_otda_jcpot_002.png
+    :class: sphx-glr-single-img
+
+
+.. rst-class:: sphx-glr-script-out
+
+ Out:
+
+ .. code-block:: none
+
+
+    (-1.85, 5.85, -4.11901398007908, 4.201462272227509)
+
+
+
+Instantiate JCPOT adaptation algorithm and fit it
+----------------------------------------------------------------------------
+
+
+.. code-block:: default
+
+    otda = ot.da.JCPOTTransport(reg_e=1, max_iter=1000, metric='sqeuclidean', tol=1e-9, verbose=True, log=True)
+    otda.fit(all_Xr, all_Yr, xt)
+
+    ws1 = otda.proportions_.dot(otda.log_['D2'][0])
+    ws2 = otda.proportions_.dot(otda.log_['D2'][1])
+
+    pl.figure(3)
+    pl.clf()
+    plot_ax(dec1, 'Source 1')
+    plot_ax(dec2, 'Source 2')
+    plot_ax(dect, 'Target')
+    print_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)
+    print_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)
+    pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+    pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+    pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+    pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+    pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+    pl.title('OT with prop estimation ({:1.3f},{:1.3f})'.format(otda.proportions_[0], otda.proportions_[1]))
+
+    pl.legend()
+    pl.axis('equal')
+    pl.axis('off')
+
+
+
+
+.. image:: /auto_examples/images/sphx_glr_plot_otda_jcpot_003.png
+    :class: sphx-glr-single-img
+
+
+.. rst-class:: sphx-glr-script-out
+
+ Out:
+
+ .. code-block:: none
+
+
+    (-1.85, 5.85, -4.11901398007908, 4.201462272227509)
+
+
+
+Run oracle transport algorithm with known proportions
+----------------------------------------------------------------------------
+
+
+.. code-block:: default
+
+    h_res = np.array([1 - pt, pt])
+
+    ws1 = h_res.dot(otda.log_['D2'][0])
+    ws2 = h_res.dot(otda.log_['D2'][1])
+
+    pl.figure(4)
+    pl.clf()
+    plot_ax(dec1, 'Source 1')
+    plot_ax(dec2, 'Source 2')
+    plot_ax(dect, 'Target')
+    print_G(ot.bregman.sinkhorn(ws1, [], otda.log_['M'][0], reg=1e-1), xs1, ys1, xt)
+    print_G(ot.bregman.sinkhorn(ws2, [], otda.log_['M'][1], reg=1e-1), xs2, ys2, xt)
+    pl.scatter(xs1[:, 0], xs1[:, 1], c=ys1, s=35, marker='x', cmap='Set1', vmax=9)
+    pl.scatter(xs2[:, 0], xs2[:, 1], c=ys2, s=35, marker='+', cmap='Set1', vmax=9)
+    pl.scatter(xt[:, 0], xt[:, 1], c=yt, s=35, marker='o', cmap='Set1', vmax=9)
+
+    pl.plot([], [], 'r', alpha=.2, label='Mass from Class 1')
+    pl.plot([], [], 'b', alpha=.2, label='Mass from Class 2')
+
+    pl.title('OT with known proportion ({:1.1f},{:1.1f})'.format(h_res[0], h_res[1]))
+
+    pl.legend()
+    pl.axis('equal')
+    pl.axis('off')
+    pl.show()
+
+
+
+.. image:: /auto_examples/images/sphx_glr_plot_otda_jcpot_004.png
+    :class: sphx-glr-single-img
+
+
+.. rst-class:: sphx-glr-script-out
+
+ Out:
+
+ .. code-block:: none
+
+    /home/rflamary/PYTHON/POT/examples/plot_otda_jcpot.py:171: UserWarning: Matplotlib is currently using agg, which is a non-GUI backend, so cannot show the figure.
+      pl.show()
+
+
+
+
+
+.. rst-class:: sphx-glr-timing
+
+   **Total running time of the script:** ( 0 minutes  4.725 seconds)
+
+
+.. _sphx_glr_download_auto_examples_plot_otda_jcpot.py:
+
+
+.. only :: html
+
+ .. container:: sphx-glr-footer
+    :class: sphx-glr-footer-example
+
+
+
+  .. container:: sphx-glr-download sphx-glr-download-python
+
+     :download:`Download Python source code: plot_otda_jcpot.py <plot_otda_jcpot.py>`
+
+
+
+  .. container:: sphx-glr-download sphx-glr-download-jupyter
+
+     :download:`Download Jupyter notebook: plot_otda_jcpot.ipynb <plot_otda_jcpot.ipynb>`
+
+
+.. only:: html
+
+ .. rst-class:: sphx-glr-signature
+
+    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_
diff --git a/docs/source/auto_examples/plot_partial_wass_and_gromov.ipynb b/docs/source/auto_examples/plot_partial_wass_and_gromov.ipynb
new file mode 100644
index 0000000..0f69ec1
--- /dev/null
+++ b/docs/source/auto_examples/plot_partial_wass_and_gromov.ipynb
@@ -0,0 +1,126 @@
+{
+  "cells": [
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "%matplotlib inline"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "\n# Partial Wasserstein and Gromov-Wasserstein example\n\n\nThis example is designed to show how to use the Partial (Gromov-)Wassertsein\ndistance computation in POT.\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "# Author: Laetitia Chapel <laetitia.chapel@irisa.fr>\n# License: MIT License\n\n# necessary for 3d plot even if not used\nfrom mpl_toolkits.mplot3d import Axes3D  # noqa\nimport scipy as sp\nimport numpy as np\nimport matplotlib.pylab as pl\nimport ot"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Sample two 2D Gaussian distributions and plot them\n--------------------------------------------------\n\nFor demonstration purpose, we sample two Gaussian distributions in 2-d\nspaces and add some random noise.\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "n_samples = 20  # nb samples (gaussian)\nn_noise = 20  # nb of samples (noise)\n\nmu = np.array([0, 0])\ncov = np.array([[1, 0], [0, 2]])\n\nxs = ot.datasets.make_2D_samples_gauss(n_samples, mu, cov)\nxs = np.append(xs, (np.random.rand(n_noise, 2) + 1) * 4).reshape((-1, 2))\nxt = ot.datasets.make_2D_samples_gauss(n_samples, mu, cov)\nxt = np.append(xt, (np.random.rand(n_noise, 2) + 1) * -3).reshape((-1, 2))\n\nM = sp.spatial.distance.cdist(xs, xt)\n\nfig = pl.figure()\nax1 = fig.add_subplot(131)\nax1.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\nax2 = fig.add_subplot(132)\nax2.scatter(xt[:, 0], xt[:, 1], color='r')\nax3 = fig.add_subplot(133)\nax3.imshow(M)\npl.show()"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Compute partial Wasserstein plans and distance,\nby transporting 50% of the mass\n----------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "p = ot.unif(n_samples + n_noise)\nq = ot.unif(n_samples + n_noise)\n\nw0, log0 = ot.partial.partial_wasserstein(p, q, M, m=0.5, log=True)\nw, log = ot.partial.entropic_partial_wasserstein(p, q, M, reg=0.1, m=0.5,\n                                                 log=True)\n\nprint('Partial Wasserstein distance (m = 0.5): ' + str(log0['partial_w_dist']))\nprint('Entropic partial Wasserstein distance (m = 0.5): ' +\n      str(log['partial_w_dist']))\n\npl.figure(1, (10, 5))\npl.subplot(1, 2, 1)\npl.imshow(w0, cmap='jet')\npl.title('Partial Wasserstein')\npl.subplot(1, 2, 2)\npl.imshow(w, cmap='jet')\npl.title('Entropic partial Wasserstein')\npl.show()"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Sample one 2D and 3D Gaussian distributions and plot them\n---------------------------------------------------------\n\nThe Gromov-Wasserstein distance allows to compute distances with samples that\ndo not belong to the same metric space. For demonstration purpose, we sample\ntwo Gaussian distributions in 2- and 3-dimensional spaces.\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "n_samples = 20  # nb samples\nn_noise = 10  # nb of samples (noise)\n\np = ot.unif(n_samples + n_noise)\nq = ot.unif(n_samples + n_noise)\n\nmu_s = np.array([0, 0])\ncov_s = np.array([[1, 0], [0, 1]])\n\nmu_t = np.array([0, 0, 0])\ncov_t = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]])\n\n\nxs = ot.datasets.make_2D_samples_gauss(n_samples, mu_s, cov_s)\nxs = np.concatenate((xs, ((np.random.rand(n_noise, 2) + 1) * 4)), axis=0)\nP = sp.linalg.sqrtm(cov_t)\nxt = np.random.randn(n_samples, 3).dot(P) + mu_t\nxt = np.concatenate((xt, ((np.random.rand(n_noise, 3) + 1) * 10)), axis=0)\n\nfig = pl.figure()\nax1 = fig.add_subplot(121)\nax1.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\nax2 = fig.add_subplot(122, projection='3d')\nax2.scatter(xt[:, 0], xt[:, 1], xt[:, 2], color='r')\npl.show()"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Compute partial Gromov-Wasserstein plans and distance,\nby transporting 100% and 2/3 of the mass\n-----------------------------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "C1 = sp.spatial.distance.cdist(xs, xs)\nC2 = sp.spatial.distance.cdist(xt, xt)\n\nprint('-----m = 1')\nm = 1\nres0, log0 = ot.partial.partial_gromov_wasserstein(C1, C2, p, q, m=m,\n                                                   log=True)\nres, log = ot.partial.entropic_partial_gromov_wasserstein(C1, C2, p, q, 10,\n                                                          m=m, log=True)\n\nprint('Partial Wasserstein distance (m = 1): ' + str(log0['partial_gw_dist']))\nprint('Entropic partial Wasserstein distance (m = 1): ' +\n      str(log['partial_gw_dist']))\n\npl.figure(1, (10, 5))\npl.title(\"mass to be transported m = 1\")\npl.subplot(1, 2, 1)\npl.imshow(res0, cmap='jet')\npl.title('Partial Wasserstein')\npl.subplot(1, 2, 2)\npl.imshow(res, cmap='jet')\npl.title('Entropic partial Wasserstein')\npl.show()\n\nprint('-----m = 2/3')\nm = 2 / 3\nres0, log0 = ot.partial.partial_gromov_wasserstein(C1, C2, p, q, m=m, log=True)\nres, log = ot.partial.entropic_partial_gromov_wasserstein(C1, C2, p, q, 10,\n                                                          m=m, log=True)\n\nprint('Partial Wasserstein distance (m = 2/3): ' +\n      str(log0['partial_gw_dist']))\nprint('Entropic partial Wasserstein distance (m = 2/3): ' +\n      str(log['partial_gw_dist']))\n\npl.figure(1, (10, 5))\npl.title(\"mass to be transported m = 2/3\")\npl.subplot(1, 2, 1)\npl.imshow(res0, cmap='jet')\npl.title('Partial Wasserstein')\npl.subplot(1, 2, 2)\npl.imshow(res, cmap='jet')\npl.title('Entropic partial Wasserstein')\npl.show()"
+      ]
+    }
+  ],
+  "metadata": {
+    "kernelspec": {
+      "display_name": "Python 3",
+      "language": "python",
+      "name": "python3"
+    },
+    "language_info": {
+      "codemirror_mode": {
+        "name": "ipython",
+        "version": 3
+      },
+      "file_extension": ".py",
+      "mimetype": "text/x-python",
+      "name": "python",
+      "nbconvert_exporter": "python",
+      "pygments_lexer": "ipython3",
+      "version": "3.6.9"
+    }
+  },
+  "nbformat": 4,
+  "nbformat_minor": 0
+}
\ No newline at end of file
diff --git a/docs/source/auto_examples/plot_partial_wass_and_gromov.py b/docs/source/auto_examples/plot_partial_wass_and_gromov.py
new file mode 100644
index 0000000..01141f2
--- /dev/null
+++ b/docs/source/auto_examples/plot_partial_wass_and_gromov.py
@@ -0,0 +1,165 @@
+# -*- coding: utf-8 -*-

+"""

+==========================

+Partial Wasserstein and Gromov-Wasserstein example

+==========================

+

+This example is designed to show how to use the Partial (Gromov-)Wassertsein

+distance computation in POT.

+"""

+

+# Author: Laetitia Chapel <laetitia.chapel@irisa.fr>

+# License: MIT License

+

+# necessary for 3d plot even if not used

+from mpl_toolkits.mplot3d import Axes3D  # noqa

+import scipy as sp

+import numpy as np

+import matplotlib.pylab as pl

+import ot

+

+

+#############################################################################

+#

+# Sample two 2D Gaussian distributions and plot them

+# --------------------------------------------------

+#

+# For demonstration purpose, we sample two Gaussian distributions in 2-d

+# spaces and add some random noise.

+

+

+n_samples = 20  # nb samples (gaussian)

+n_noise = 20  # nb of samples (noise)

+

+mu = np.array([0, 0])

+cov = np.array([[1, 0], [0, 2]])

+

+xs = ot.datasets.make_2D_samples_gauss(n_samples, mu, cov)

+xs = np.append(xs, (np.random.rand(n_noise, 2) + 1) * 4).reshape((-1, 2))

+xt = ot.datasets.make_2D_samples_gauss(n_samples, mu, cov)

+xt = np.append(xt, (np.random.rand(n_noise, 2) + 1) * -3).reshape((-1, 2))

+

+M = sp.spatial.distance.cdist(xs, xt)

+

+fig = pl.figure()

+ax1 = fig.add_subplot(131)

+ax1.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')

+ax2 = fig.add_subplot(132)

+ax2.scatter(xt[:, 0], xt[:, 1], color='r')

+ax3 = fig.add_subplot(133)

+ax3.imshow(M)

+pl.show()

+

+#############################################################################

+#

+# Compute partial Wasserstein plans and distance,

+# by transporting 50% of the mass

+# ----------------------------------------------

+

+p = ot.unif(n_samples + n_noise)

+q = ot.unif(n_samples + n_noise)

+

+w0, log0 = ot.partial.partial_wasserstein(p, q, M, m=0.5, log=True)

+w, log = ot.partial.entropic_partial_wasserstein(p, q, M, reg=0.1, m=0.5,

+                                                 log=True)

+

+print('Partial Wasserstein distance (m = 0.5): ' + str(log0['partial_w_dist']))

+print('Entropic partial Wasserstein distance (m = 0.5): ' +

+      str(log['partial_w_dist']))

+

+pl.figure(1, (10, 5))

+pl.subplot(1, 2, 1)

+pl.imshow(w0, cmap='jet')

+pl.title('Partial Wasserstein')

+pl.subplot(1, 2, 2)

+pl.imshow(w, cmap='jet')

+pl.title('Entropic partial Wasserstein')

+pl.show()

+

+

+#############################################################################

+#

+# Sample one 2D and 3D Gaussian distributions and plot them

+# ---------------------------------------------------------

+#

+# The Gromov-Wasserstein distance allows to compute distances with samples that

+# do not belong to the same metric space. For demonstration purpose, we sample

+# two Gaussian distributions in 2- and 3-dimensional spaces.

+

+n_samples = 20  # nb samples

+n_noise = 10  # nb of samples (noise)

+

+p = ot.unif(n_samples + n_noise)

+q = ot.unif(n_samples + n_noise)

+

+mu_s = np.array([0, 0])

+cov_s = np.array([[1, 0], [0, 1]])

+

+mu_t = np.array([0, 0, 0])

+cov_t = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]])

+

+

+xs = ot.datasets.make_2D_samples_gauss(n_samples, mu_s, cov_s)

+xs = np.concatenate((xs, ((np.random.rand(n_noise, 2) + 1) * 4)), axis=0)

+P = sp.linalg.sqrtm(cov_t)

+xt = np.random.randn(n_samples, 3).dot(P) + mu_t

+xt = np.concatenate((xt, ((np.random.rand(n_noise, 3) + 1) * 10)), axis=0)

+

+fig = pl.figure()

+ax1 = fig.add_subplot(121)

+ax1.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')

+ax2 = fig.add_subplot(122, projection='3d')

+ax2.scatter(xt[:, 0], xt[:, 1], xt[:, 2], color='r')

+pl.show()

+

+

+#############################################################################

+#

+# Compute partial Gromov-Wasserstein plans and distance,

+# by transporting 100% and 2/3 of the mass

+# -----------------------------------------------------

+

+C1 = sp.spatial.distance.cdist(xs, xs)

+C2 = sp.spatial.distance.cdist(xt, xt)

+

+print('-----m = 1')

+m = 1

+res0, log0 = ot.partial.partial_gromov_wasserstein(C1, C2, p, q, m=m,

+                                                   log=True)

+res, log = ot.partial.entropic_partial_gromov_wasserstein(C1, C2, p, q, 10,

+                                                          m=m, log=True)

+

+print('Partial Wasserstein distance (m = 1): ' + str(log0['partial_gw_dist']))

+print('Entropic partial Wasserstein distance (m = 1): ' +

+      str(log['partial_gw_dist']))

+

+pl.figure(1, (10, 5))

+pl.title("mass to be transported m = 1")

+pl.subplot(1, 2, 1)

+pl.imshow(res0, cmap='jet')

+pl.title('Partial Wasserstein')

+pl.subplot(1, 2, 2)

+pl.imshow(res, cmap='jet')

+pl.title('Entropic partial Wasserstein')

+pl.show()

+

+print('-----m = 2/3')

+m = 2 / 3

+res0, log0 = ot.partial.partial_gromov_wasserstein(C1, C2, p, q, m=m, log=True)

+res, log = ot.partial.entropic_partial_gromov_wasserstein(C1, C2, p, q, 10,

+                                                          m=m, log=True)

+

+print('Partial Wasserstein distance (m = 2/3): ' +

+      str(log0['partial_gw_dist']))

+print('Entropic partial Wasserstein distance (m = 2/3): ' +

+      str(log['partial_gw_dist']))

+

+pl.figure(1, (10, 5))

+pl.title("mass to be transported m = 2/3")

+pl.subplot(1, 2, 1)

+pl.imshow(res0, cmap='jet')

+pl.title('Partial Wasserstein')

+pl.subplot(1, 2, 2)

+pl.imshow(res, cmap='jet')

+pl.title('Entropic partial Wasserstein')

+pl.show()

diff --git a/docs/source/auto_examples/plot_screenkhorn_1D.ipynb b/docs/source/auto_examples/plot_screenkhorn_1D.ipynb
new file mode 100644
index 0000000..1c27d3b
--- /dev/null
+++ b/docs/source/auto_examples/plot_screenkhorn_1D.ipynb
@@ -0,0 +1,108 @@
+{
+  "cells": [
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "%matplotlib inline"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "\n# 1D Screened optimal transport\n\n\nThis example illustrates the computation of Screenkhorn:\nScreening Sinkhorn Algorithm for Optimal transport.\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "# Author: Mokhtar Z. Alaya <mokhtarzahdi.alaya@gmail.com>\n#\n# License: MIT License\n\nimport numpy as np\nimport matplotlib.pylab as pl\nimport ot.plot\nfrom ot.datasets import make_1D_gauss as gauss\nfrom ot.bregman import screenkhorn"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Generate data\n-------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "n = 100  # nb bins\n\n# bin positions\nx = np.arange(n, dtype=np.float64)\n\n# Gaussian distributions\na = gauss(n, m=20, s=5)  # m= mean, s= std\nb = gauss(n, m=60, s=10)\n\n# loss matrix\nM = ot.dist(x.reshape((n, 1)), x.reshape((n, 1)))\nM /= M.max()"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Plot distributions and loss matrix\n----------------------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "pl.figure(1, figsize=(6.4, 3))\npl.plot(x, a, 'b', label='Source distribution')\npl.plot(x, b, 'r', label='Target distribution')\npl.legend()\n\n# plot distributions and loss matrix\n\npl.figure(2, figsize=(5, 5))\not.plot.plot1D_mat(a, b, M, 'Cost matrix M')"
+      ]
+    },
+    {
+      "cell_type": "markdown",
+      "metadata": {},
+      "source": [
+        "Solve Screenkhorn\n-----------------------\n\n"
+      ]
+    },
+    {
+      "cell_type": "code",
+      "execution_count": null,
+      "metadata": {
+        "collapsed": false
+      },
+      "outputs": [],
+      "source": [
+        "# Screenkhorn\nlambd = 2e-03  # entropy parameter\nns_budget = 30  # budget number of points to be keeped in the source distribution\nnt_budget = 30  # budget number of points to be keeped in the target distribution\n\nG_screen = screenkhorn(a, b, M, lambd, ns_budget, nt_budget, uniform=False, restricted=True, verbose=True)\npl.figure(4, figsize=(5, 5))\not.plot.plot1D_mat(a, b, G_screen, 'OT matrix Screenkhorn')\npl.show()"
+      ]
+    }
+  ],
+  "metadata": {
+    "kernelspec": {
+      "display_name": "Python 3",
+      "language": "python",
+      "name": "python3"
+    },
+    "language_info": {
+      "codemirror_mode": {
+        "name": "ipython",
+        "version": 3
+      },
+      "file_extension": ".py",
+      "mimetype": "text/x-python",
+      "name": "python",
+      "nbconvert_exporter": "python",
+      "pygments_lexer": "ipython3",
+      "version": "3.6.9"
+    }
+  },
+  "nbformat": 4,
+  "nbformat_minor": 0
+}
\ No newline at end of file
diff --git a/docs/source/auto_examples/plot_screenkhorn_1D.py b/docs/source/auto_examples/plot_screenkhorn_1D.py
new file mode 100644
index 0000000..840ead8
--- /dev/null
+++ b/docs/source/auto_examples/plot_screenkhorn_1D.py
@@ -0,0 +1,68 @@
+# -*- coding: utf-8 -*-
+"""
+===============================
+1D Screened optimal transport
+===============================
+
+This example illustrates the computation of Screenkhorn:
+Screening Sinkhorn Algorithm for Optimal transport.
+"""
+
+# Author: Mokhtar Z. Alaya <mokhtarzahdi.alaya@gmail.com>
+#
+# License: MIT License
+
+import numpy as np
+import matplotlib.pylab as pl
+import ot.plot
+from ot.datasets import make_1D_gauss as gauss
+from ot.bregman import screenkhorn
+
+##############################################################################
+# Generate data
+# -------------
+
+#%% parameters
+
+n = 100  # nb bins
+
+# bin positions
+x = np.arange(n, dtype=np.float64)
+
+# Gaussian distributions
+a = gauss(n, m=20, s=5)  # m= mean, s= std
+b = gauss(n, m=60, s=10)
+
+# loss matrix
+M = ot.dist(x.reshape((n, 1)), x.reshape((n, 1)))
+M /= M.max()
+
+##############################################################################
+# Plot distributions and loss matrix
+# ----------------------------------
+
+#%% plot the distributions
+
+pl.figure(1, figsize=(6.4, 3))
+pl.plot(x, a, 'b', label='Source distribution')
+pl.plot(x, b, 'r', label='Target distribution')
+pl.legend()
+
+# plot distributions and loss matrix
+
+pl.figure(2, figsize=(5, 5))
+ot.plot.plot1D_mat(a, b, M, 'Cost matrix M')
+
+##############################################################################
+# Solve Screenkhorn
+# -----------------------
+
+# Screenkhorn
+lambd = 2e-03  # entropy parameter
+ns_budget = 30  # budget number of points to be keeped in the source distribution
+nt_budget = 30  # budget number of points to be keeped in the target distribution
+
+G_screen = screenkhorn(a, b, M, lambd, ns_budget, nt_budget, uniform=False, restricted=True, verbose=True)
+pl.figure(4, figsize=(5, 5))
+ot.plot.plot1D_mat(a, b, G_screen, 'OT matrix Screenkhorn')
+pl.show()