Week 3. How Moving Earth Changed GANs!.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Recap: $f$-divergences:\n",
    "In the [last note](Week%202.%20A%20General%20Framework%20for%20GANs.ipynb) we saw how [[1]](#References) generalized Generative Adversarial Networks (GANs) to approximately minimize any $f$-divergence. Concretely, theorem 1 of [1] gave assumptions for which we can approximately minimize \n",
    "\n",
    "$$D_f(p_d || p_g) = \\int_x p_g(x) f\\left(\\frac{p_d(x)}{p_g(x)}\\right) dx$$\n",
    "\n",
    "by playing the following GAN game\n",
    "\n",
    "$$\n",
    "\\min_{G_\\theta} \\max_{D_\\omega} F(\\omega, \\theta)\n",
    "\\quad\\text{where}\\quad F(\\omega, \\theta) = \\mathbb{E}_{x\\sim p_d}[{D_\\omega}(x)] - \\mathbb{E}_{z \\sim p_z} [f^\\ast({D_\\omega}({G_\\theta}(z)))]\n",
    "$$\n",
    "\n",
    "Here $f^*$ is the complex conjungate. Unfortunately, many $f$-divergences become very large for non-overlapping distributions. In some cases they are not even continous with respect to the model parameter (and thus not differentiable). Since we train all our models with gradient descent this is quite troubling. \n",
    "\n",
    "<hr>\n",
    "<b>Example 1: </b> (Uniform Distribution, see [2])<br>\n",
    "Let $Z\\sim U(0,1)$ and $A=(0, Z)$ and $B=(\\theta, Z)$. The resulting distributions are visualized below for $\\Theta=0, 0.5, ..., 5$: \n",
    "\n",
    "<video autoplay src=\"gifs/uniform_example/animation.mp4\" controls loop style=\"width: 80%; margin-left:10%; margin-top: 12px;\">\n",
    "Your browser does not support the video tag\n",
    "</video>\n",
    "<p style=\"width: 80%; margin-left: 10%; text-align: center;\"><b>Figure 1:</b> Uniform distributions visualized.</p>\n",
    "\n",
    "For a fixed $\\theta$ we can think of $A$ as our initial attempt at learning $B$. The error of $A$ is then captured by our $f$-divergence, for simplicity we consider only the KL divergence:  \n",
    "\n",
    "$$ D_{kl}(A, B)= \\int_x A(x) f\\left(\\frac{A(x)}{B(x)}\\right) dx$$\n",
    "\n",
    "When $\\theta=0$ the distributions are equal $A=B$ and $D_{kl}(A, A)=0$. But when $\\theta\\neq 0$ the $x$'s where $A(x)\\neq 0$ we have that $B(x)=0$. The KL divergence then becomes: \n",
    "\n",
    "$$ D_{kl}(A, B)=\\begin{cases}0 &\\text{if } \\theta=0\\\\\\infty & \\text{else}\\end{cases}$$\n",
    "\n",
    "This is not continous in $\\theta$ (and thus not differentiable). In [[3]](#References) they introduce a new distance measure $W(A,B)$ for which it holds that $W(A,B)=|\\theta|$ which is continuous and differentiable almost everywhere. \n",
    "\n",
    "<hr>\n",
    "\n",
    "The above example highlights a serious problem with $f$-divergence for non-overlapping distributions. That said, the example is perfectly non-overlapping. Is the problem also present for normal distributions with different means? Normal distributions are defined on all of $\\mathbb{R}$ so they will always overlap, even though the \"size\" of the overlap decreases rapidly as the means move further apart. \n",
    "\n",
    "<hr>\n",
    "\n",
    "<b>Example 2: </b>(Normal Distributions)<br>\n",
    "The KL divergence of normal distributions have a closed form given below: \n",
    "\n",
    "<b>Lemma: </b>Let $p=N(\\mu_1,\\sigma_1)$ and $q=N(\\mu_2, \\sigma_2)$ then $D_{kl}(p||q) = \\log \\frac{\\sigma_2}{\\sigma_1}+\\frac{\\sigma_1^2+(\\mu_1-\\mu_2)^2}{2\\sigma_2^2}-\\frac{1}{2}$. <br>\n",
    "<b>Proof: </b> See [[3]](#References). \n",
    "\n",
    "<b>Corollary: </b>If $A=N(0,1)$ and $B=N(\\mu, 1)$ then $D_{kl}(A||B)=\\frac{\\mu^2}{2}$.  \n",
    "<b>Proof: </b> Substitute $\\sigma_1=\\sigma_2=1$ and $\\mu_1=1$ in the lemma above. \n",
    "\n",
    "The code below visualizes how the KL divergence of $N(0,1)$ and $N(\\mu, 1)$ increase with $\\mu$ by using the above Corollary. \n",
    "\n",
    "<video src=\"gifs/normal_kl/animation.mp4\" controls autoplay loop style=\"width: 80%; margin-left:10%; margin-top: 12px;\">\n",
    "Your browser does not support the video tag\n",
    "</video>\n",
    "<p style=\"width: 80%; margin-left: 10%; text-align: center;\"><b>Figure 2:</b> Visualization of KL divergence of $N(0,1)$ and $N(\\mu, 1)$.</p>\n",
    "\n",
    "In this case the KL divergence is continous and differentiable in $\\mu$. If we consider a function $g(x)=ax+b$ that attempts to transform $A$ into $B$ the KL divergence of $g(A)$ and $B$ would be both continuous and differentiable with respect to the model parameters (a,b). \n",
    "\n",
    "In other words: for normal distributions the KL divergence is somehow well behaved, even though it increases quadratically with $\\mu$. The following post on stackexchange prove an upperbound and lowerbound on the distance measure $W(A, B)$ for normal distributions. The upperbound is asymptotically the same as the KL divergence but the lowerbound is asymptotically better. We leave it as an exercise to make a tight bound (this should only be attempted to solve after reading the entire notebook). Is the $W$ distance smaller than KL in terms of $\\mu$?\n",
    "\n",
    "https://stats.stackexchange.com/questions/83741/earth-movers-distance-emd-between-two-gaussians\n",
    "\n",
    "<hr>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Moving Earth\n",
    "Consider the following generated distribution $p_g$ which we want to turn into $p_d$: \n",
    "\n",
    "<img src=\"figs/earth_mover_example.png\" />\n",
    "\n",
    "We need to define some measure of distance between these two distributions. Since we usually define $p_g$ by transforming some sorce of randomnes $p_z$ by using $G$, it might be meaningful to define our measure as \"how much we need to change $p_g$ to obtain $p_d$\". Intuitively, you can think of $p_g$ and $p_d$ above as piles of earth. One possible distance measure between these two \"piles of earth\" is then how much work we need to perform to turn the first pile of earth into the other pile. \n",
    "\n",
    "Let $\\gamma(x,y)$ denote a plan for how we want to move the earth. For example, $\\gamma(1,4)=3$ means we move $3$ units of earth from $p_g$ at $x=1$ to $x=4$. If we also let $\\gamma(2,3)=1$ our plan would have succesfully transformed $p_d$ into $p_g$, see the animation below: \n",
    "\n",
    "<video src=\"gifs/transport_plan/animation.mp4\" controls autoplay loop style=\"width: 80%; margin-left:10%; margin-top: 12px;\">\n",
    "Your browser does not support the video tag\n",
    "</video>\n",
    "<p style=\"width: 80%; margin-left: 10%; text-align: center;\"><b>Figure 3:</b> Exampe of how \"moving earth\" constitutes a measure of difference between two distributions (In this case, $p_g$ and $p_d$).</p>\n",
    "\n",
    "Furthermore, we add a cost $c(x,y)$ that tells us how much it costs to transport earth from $x$ to $y$; for example we could have $c(x,y)=||x-y||_2^2$. The total cost for a plan $\\gamma$ is then given\n",
    "\n",
    "$$\\sum_{x,y}\\gamma(x,y) c(x,y)$$\n",
    "\n",
    "Since there are many plans, we choose the distance to be the optimal plan, namely, the one that transports the earth most cost effectively\n",
    "\n",
    "$$EM(p_d, p_g)=\\min_\\gamma \\; \\sum_{x,y} \\gamma(x,y)c(x,y)$$\n",
    "\n",
    "But we are only interested in plans $\\gamma$ that transform $p_g$ into $p_d$. To satisfy this constraint, we require that all plans $\\gamma$ must satisfy the following: \n",
    "\n",
    "$$ \\sum_{m} \\gamma(m, x) = p_d(x) \\quad\\quad \\sum_{m} \\gamma(x, m)=p_g(x)$$\n",
    "\n",
    "The first requirement ensures that all earth moved to $x$ sums up to the amount $p_d$ has at position $x$, in other words, it ensures that the plan actually turns $p_g$ into $p_d$. The second requirement ensures that everything moved away from $p_g$ at position $x$ is exactly how much $p_g$ has a position $x$ (notice that $\\gamma(1,1)=1$ in our above example). \n",
    "\n",
    "The above constraints mean that $\\gamma$ should marganilize into $p_d$ or $p_g$ depending on which variable we marginalize over. If we let $P[p_d, p_g]$ denote the set of joint probability distributions that marginalizes to $p_d$ and $p_g$, we get that \n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "EM(p_d, p_g)&=\\min_{\\gamma \\in P[p_d, p_g]}\\; \\sum_{x,y} \\gamma(x,y)c(x,y)\\\\\n",
    "&=\\min_{\\gamma \\in P[p_d, p_g]}\\; \\mathbb{E}_{(x,y)\\sim \\gamma}\\left[c(x,y) \\right]\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "If we let $c(x,y)=||x-y||$ and consider continous distributions $p_d$ and $p_g$ we obtain the Wasserstein distance (we use Earth Mover and Wasserstein distance interchangeably)\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "W(p_d, p_g)&=\\inf_{\\gamma \\in P[p_d, p_g]}\\; \\mathbb{E}_{(x,y)\\sim \\gamma}\\left[||x-y||\\right]\n",
    "\\end{align}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Explaining Theorem 1 and Corollary 2\n",
    "We explain libschitz later in the end of this section: \n",
    "\n",
    "<b>Theorem 1. </b>Let $x_1,...,x_n\\sim p_d$ and $z\\sim p_z$ so $p_g$ denotes implicit distribution defined by $g_\\Theta(z)$ where $\\Theta$ are the weights of $g$. Then\n",
    "1. If $g$ is continous in $\\Theta$ then $W(p_d, p_g)$ is continous in $\\Theta$. \n",
    "2. If $g$ is also Lipschitz then $W(p_d, p_g)$ is also differentiable almost everywhere. \n",
    "3. The above statements are not true for the JS and KL divergence. \n",
    "\n",
    "In Machine Learning we often frame our problems as optimization problems and minimize error/loss by some variant of gradient descent. For this to be succesfull, we'd like our error/loss to be differentiable (and thus continuous). Sometimes (as with relu non-linearity) we settle for piecewise differentiable. \n",
    "\n",
    "In example 1 above, we saw that the KL divergence is undefined in all but one point, so this is not differentiable, it is not even continous. The takeway from Theorem 1 is that the Wasserstein GAN (under fitting assumptions) is both continous in the model parameretes and differentiable \"<a href=\"https://en.wikipedia.org/wiki/Almost_everywhere\">almost everywhere</a>\", which JS and KL is not. This provides a strong theoretical justification for using the EM distance instead of JS or KL. \n",
    "\n",
    "<b>Open Problem: </b>Does all $f$-divergences have this problem? Can you find an $f$-divergence that does not have this problem?  \n",
    "\n",
    "Since the generator is implemented as a Neural Network, one might ask if Neural Networks satisfy the wanted properties? The following corollary gives a positive answer to this question: \n",
    "\n",
    "<b>Corollary 1. </b> Feed Forward Neural Networks with $p_z$ that has bounded expectation satisfies asummption 1+2 and therefore Wasserstein is continuous everywhere and differentiable almost everywhere.  \n",
    "\n",
    "\n",
    "### When is a function $L$-libschitz?\n",
    "\n",
    "<b>Definition: </b> A function $f:\\mathcal{X}\\rightarrow \\mathbb{R}$ is $L$-libschitz if $$\\forall x,y\\in \\mathcal{X}: \\quad \\frac{||f(x)-f(y)||}{||x-y||}\\le L$$\n",
    "\n",
    "Consider for example a linear function $f(x)=Wx$ where $W$ is invertible and thus has all non-zero eigenvalues. We then get that \n",
    "\n",
    "$$||f(x)-f(y)||/||x-y|| = ||W(x-y)||/||x-y||$$\n",
    "\n",
    "A linear transform can not increase the norm of a vector by more than its maximal eigenvalue. Since $(x-y)$ is a vector, we get that $||W(x-y)||\\le \\lambda_{max}(W)||x-y||$ and thus see that\n",
    "\n",
    "$$||W(x-y)||/||x-y|| \\le \\lambda_{max}$$\n",
    "\n",
    "It thus follows that linear transforms that have all non-zero eigenvalues has the libschitz constant $\\lambda_{max}(W)$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# How do we actually minimize Earth Mover? Seems intractable. \n",
    "Ok. So we introduce the earth mover distance: \n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "W(p_d, p_g)&=\\inf_{\\gamma \\in P[p_d, p_g]}\\; \\mathbb{E}_{(x,y)\\sim \\gamma}\\left[||x-y||\\right]\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "Furthermore, we saw a theorem that showed that the earth mover distance has nicer properties than JS and KL. But the above distance looks intractable, how could we compute this function? Fortunately, Kantorovich and Rubinstein manages to show that the Wasserstein distance can be computed as\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "W(p_d, p_g)&=\\sup_{||f||_L \\le 1}\\; \\mathbb{E}_{x\\sim p_d}\\left[f(x)\\right] - \\mathbb{E}_{x\\sim p_g}\\left[f(x)\\right]\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "where $||f||_L\\le 1$ means that $f$ has a libschitz constant $1$. Since we only care about minimizing $W$ (actually we only care about gradients of $W$), changing the constraint to $||f|_L\\le k$ doesn't change much, the gradient will be scaled by a factor $k$. Thus, if we consider a limited family of function $F$ that are all $K$-lipschitz for some $K$ minimizing over this set of functions minimize the above equation up to a constant factor. \n",
    "\n",
    "We would like our set of functions $F$ to be specified by a certain Neural Network architecture. To make sure this family of neural networks has a bounded libschitz constant, we could just clip the weights of each layer. Intuitively, if we consider a single layer, its maximum eigenvalue must be bounded; if one accepts that this holds under composition with non-linearities (highly nontrivial), we would get that weight clipping constrains the family to have a bounded libschitz constant. \n",
    "\n",
    "A \"theoretically cleaner\" approach was later introduced by [[4]](#References) who suggest spectral normalization. This technique is the topic of week 4. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Code for moving normal distributions\n",
    "The following two peaces of code is what generated Figure 1 and Figure 2, respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt \n",
    "from scipy.stats import norm\n",
    "\n",
    "plt.ion()\n",
    "fig, ax = plt.subplots(1,figsize=(10, 2))\n",
    "fig.suptitle('Uniform and Moving Uniform Distribution', fontsize=12)\n",
    "\n",
    "\n",
    "n = 25\n",
    "theta=5\n",
    "\n",
    "for i in range(11): \n",
    "    theta = i / 2\n",
    "    p = np.concatenate( (np.zeros((n,1 )), np.random.uniform(0, 1, size=(n,1))), axis=1)\n",
    "    q = np.concatenate( (theta*np.ones((n,1 )), np.random.uniform(0, 1, size=(n,1))), axis=1)\n",
    "\n",
    "    ax.cla()\n",
    "    ax.plot(p[:, 0], p[:, 1], 'x', label=\"A=(0, Uniform)\")\n",
    "    ax.plot(q[:, 0], q[:, 1], 'o', label=\"B=(%.1f, Uniform)\"%theta)\n",
    "    ax.set_ylim([0, 1])\n",
    "    ax.set_xlim([-0.5, 5.5])\n",
    "    ax.legend(loc=1)\n",
    "\n",
    "    plt.savefig(\"gifs/uniform_example/%i.png\"%i)\n",
    "    fig.canvas.draw()\n",
    "    plt.pause(.1)\n",
    "    \n",
    "    \n",
    "for i in range(11): \n",
    "    theta = 5 - i / 2\n",
    "    p = np.concatenate( (np.zeros((n,1 )), np.random.uniform(0, 1, size=(n,1))), axis=1)\n",
    "    q = np.concatenate( (theta*np.ones((n,1 )), np.random.uniform(0, 1, size=(n,1))), axis=1)\n",
    "\n",
    "    ax.cla()\n",
    "    ax.plot(p[:, 0], p[:, 1], 'x', label=\"A=(0, Uniform)\")\n",
    "    ax.plot(q[:, 0], q[:, 1], 'o', label=\"B=(%.1f, Uniform)\"%theta)\n",
    "    ax.set_ylim([0, 1])\n",
    "    ax.set_xlim([-0.5, 5.5])\n",
    "    ax.legend(loc=1)\n",
    "\n",
    "\n",
    "    plt.savefig(\"gifs/uniform_example/%i.png\"%(i+10))\n",
    "    fig.canvas.draw()\n",
    "    plt.pause(.1)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt \n",
    "from scipy.stats import norm\n",
    "\n",
    "def kl_div(mu1, std1, mu2, std2): \n",
    "    #..https://stats.stackexchange.com/questions/7440/kl-divergence-between-two-univariate-gaussians\n",
    "    return np.log(std1 / std2) + (std1**2 + (mu1 - mu2)**2) / (2*std2**2) - 1/2\n",
    "\n",
    "# compute kl divergence between N(0, 1) and N(mu, 1) as mu increase. \n",
    "xs = np.arange(0,  30, 1)\n",
    "xs_ = np.arange(-5, 30, 0.1)\n",
    "ys = kl_div(0, 1, xs, 1)\n",
    "\n",
    "plt.ion()\n",
    "fig, ax = plt.subplots(1, 2, figsize=(10,3))\n",
    "#ax.set_xlabel(\"increasing mean\")\n",
    "fig.suptitle('KL of N(0, 1) and Moving Normal Distribution', fontsize=12)\n",
    "\n",
    "# plot it moving from left to right. \n",
    "for i in range(0, 15):     \n",
    "    ax[0].cla()\n",
    "    ax[1].cla()\n",
    "    ax[0].plot(xs_, norm.pdf(xs_, 0, 1), label=\"N(0, 1)\", alpha=0.3, lw=4)\n",
    "    ax[0].plot(xs_, norm.pdf(xs_, 0+i, 1), label=\"N(%i, 1)\"%int(i), alpha=0.3)\n",
    "    ax[1].plot(xs[:i], ys[:i], label=\"KL\")\n",
    "    ax[0].plot([-5, 20], [ys[i], ys[i]], label=\"KL\")\n",
    "    ax[0].set_yscale('log')\n",
    "    ax[0].set_ylim([0.1, 1000])\n",
    "    ax[1].set_ylim([0, 100])\n",
    "    ax[1].set_xlim([0, 20])\n",
    "    ax[0].set_xlim([-5, 20])\n",
    "    ax[0].legend(loc=1)\n",
    "    ax[1].legend(loc=1)\n",
    "    plt.savefig(\"gifs/normal_kl/%i.png\"%i)\n",
    "    fig.canvas.draw()\n",
    "    plt.pause(.1)\n",
    "    \n",
    "for i in range(15, -1, -1):\n",
    "    ax[0].cla()\n",
    "    ax[1].cla()\n",
    "    ax[0].plot(xs_, norm.pdf(xs_, 0, 1), label=\"N(0, 1)\", alpha=0.3, lw=4)\n",
    "    ax[0].plot(xs_, norm.pdf(xs_, 0+i, 1), label=\"N(%i, 1)\"%int(i), alpha=0.3)\n",
    "    ax[1].plot(xs[:i], ys[:i], label=\"KL\")\n",
    "    ax[0].plot([-5, 20], [ys[i], ys[i]], label=\"KL\")\n",
    "    ax[0].set_yscale('log')\n",
    "    ax[0].set_ylim([0.1, 1000])\n",
    "    ax[1].set_ylim([0, 100])\n",
    "    ax[1].set_xlim([0, 20])\n",
    "    ax[0].set_xlim([-5, 20])\n",
    "\n",
    "    ax[0].legend(loc=1)\n",
    "    ax[1].legend(loc=1)\n",
    "    plt.savefig(\"gifs/normal_kl/%i.png\"%(30-i))\n",
    "    fig.canvas.draw()\n",
    "    plt.pause(.1)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The code for generating Figure 3 can be found below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "p_g = np.array([1,1,1,1,2,2,2,3,3,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "fig, ax = plt.subplots(1, 2, figsize=(10,2))\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "plt.savefig(\"figs/earth_mover_example.png\")\n",
    "plt.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.ion()\n",
    "fig, ax = plt.subplots(1, 2, figsize=(10,2))\n",
    "\n",
    "p_g = np.array([1,1,1,1,2,2,2,3,3,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "ax[0].cla()\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "ax[0].set_ylim([0,4])\n",
    "plt.savefig(\"gifs/transport_plan/0.png\")\n",
    "fig.canvas.draw()\n",
    "plt.pause(1)\n",
    "\n",
    "p_g = np.array([1, 2,2,2,3,3,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "ax[0].cla()\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "ax[0].set_ylim([0,4])\n",
    "plt.savefig(\"gifs/transport_plan/1.png\")\n",
    "fig.canvas.draw()\n",
    "plt.pause(1)\n",
    "\n",
    "p_g = np.array([1,2,2,2,3,3,4,4,4,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "ax[0].cla()\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "ax[0].set_ylim([0,4])\n",
    "plt.savefig(\"gifs/transport_plan/2.png\")\n",
    "fig.canvas.draw()\n",
    "plt.pause(1)\n",
    "\n",
    "p_g = np.array([1,2,2,3,3,4,4,4,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "ax[0].cla()\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "ax[0].set_ylim([0,4])\n",
    "plt.savefig(\"gifs/transport_plan/3.png\")\n",
    "fig.canvas.draw()\n",
    "plt.pause(1)\n",
    "\n",
    "p_g = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "p_d = np.array([1,2,2,3,3,3,4,4,4,4])\n",
    "\n",
    "ax[0].cla()\n",
    "ax[0].hist(p_g, color=\"C0\")\n",
    "ax[0].set_title(\"$p_g$\")\n",
    "ax[1].hist(p_d, color=\"C1\")\n",
    "ax[1].set_title(\"$p_d$\")\n",
    "ax[0].set_ylim([0,4])\n",
    "plt.savefig(\"gifs/transport_plan/4.png\")\n",
    "fig.canvas.draw()\n",
    "plt.pause(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# References \n",
    "[[1]](http://papers.nips.cc/paper/6066-f-gan-training-generative-neural-samplers-using-variational-divergence-minimization.pdf) Nowozin, S., Cseke, B. and Tomioka, R., 2016. f-gan: Training generative neural samplers using variational divergence minimization. In Advances in Neural Information Processing Systems (pp. 271-279).\n",
    "\n",
    "[[2]](https://arxiv.org/abs/1701.07875) Arjovsky, M., Chintala, S. and Bottou, L., 2017, July. Wasserstein generative adversarial networks. In International Conference on Machine Learning (pp. 214-223).\n",
    "\n",
    "[3] https://stats.stackexchange.com/questions/7440/kl-divergence-between-two-univariate-gaussians\n",
    "\n",
    "[[4]](https://arxiv.org/pdf/1802.05957.pdf) Miyato, T., Kataoka, T., Koyama, M. and Yoshida, Y., 2018. Spectral normalization for generative adversarial networks. arXiv preprint arXiv:1802.05957."
   ]
  }
 ],
 "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.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}