notebooks/slap_two_layer_neural_network.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# slap_two_layer_neural_network\n",
    "\n",
    "We show implementation of a two-layer neural network by [Size Linear Algebra Library (SLAP)](http://akabe.github.io/slap/).\n",
    "It is a high-level wrapper of [Lacaml](https://github.com/mmottl/lacaml), a binding of linear algebra libraries LAPACK and BLAS,\n",
    "with type-level size checking for vector and matrix operations.\n",
    "SLAP helps your debug by detecting inconsistency of dimensions\n",
    "\n",
    "- **at compile time** (instead of runtime) and\n",
    "- **at higher level** (i.e., at a caller site rather than somewhere deep inside of a call stack).\n",
    "\n",
    "For example, addition of vectors of different sizes causes a type error at compile time, and dynamic errors such as exceptions do not happen. For most high-level matrix operations, the consistency of sizes is verified statically. (Certain low-level operations, like accesses to elements by indices, need dynamic checks.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Module Archimedes loaded and aliased as A.\n"
     ]
    }
   ],
   "source": [
    "#thread ;;\n",
    "#require \"core,slap,slap.top,slap.ppx,jupyter-archimedes\" ;;\n",
    "#print_depth 3 ;; (* Avoid too deep pretty printing *)\n",
    "#print_length 30 ;; (* Avoid too long pretty printing *)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "open Core ;;\n",
    "open Slap.Common ;;\n",
    "open Slap.D ;;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dimensions, vectors, and matrics\n",
    "\n",
    "The dimensions of vectors and matrics have type `'n Slap.Size.t` like `int`,\n",
    "but type parameter `'n` is the size of itselt as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val input_dim : Slap.Size.two Slap.Size.t = 2\n"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "val hidden_dim : Slap.Size.five Slap.Size.t = 5\n"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "val output_dim : Slap.Size.one Slap.Size.t = 1\n"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let input_dim = Slap.Size.two ;;\n",
    "let hidden_dim = Slap.Size.five ;;\n",
    "let output_dim = Slap.Size.one ;;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similarly, `('n, 'a) vec` is the type of `'n`-dimensional vectors like:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val x : (Slap.Size.two, 'a) Slap_D.vec = R1 R2\n",
       "                                          1  1\n"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let x = Vec.make1 input_dim"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The second type parameter is a _contiguous or discrete_ flag. For the time being, you do not need to consider it because we only use contiguous matrices.\n",
    "\n",
    "The matrix type `('m, 'n, 'a) mat` has two dimensional type parameters: `'m`, `'n` are the numbers of rows and columns, respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val a : (Slap.Size.five, Slap.Size.two, 'a) Slap_D.mat =\n",
       "            C1        C2\n",
       "  R1  0.683384  0.683698\n",
       "  R2  0.903266 -0.321421\n",
       "  R3  0.600824  0.400497\n",
       "  R4 -0.845824  0.709542\n",
       "  R5  0.153864  0.610164\n"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let a = Mat.random hidden_dim input_dim ;;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Prediction\n",
    "\n",
    "### Single-layer neural networks\n",
    "\n",
    "Let $\\mathbf{x}$ be an input vector of a layer, $\\mathbf{W}$ be a weight matrix, and $\\mathbf{b}$ be a bias.\n",
    "Output vector $\\mathbf{y}$ is given by\n",
    "\n",
    "$$\\mathbf{y} = \\sigma\\left( \\mathbf{W} \\mathbf{x} + \\mathbf{b} \\right)$$\n",
    "\n",
    "where $\\sigma$ is a sigmoid function defined as\n",
    "\n",
    "$$\\sigma(a) = \\frac{1}{1 + \\exp(-a)}.$$\n",
    "\n",
    "A Level-2 BLAS function `gemv ~trans:normal ~alpha a x ~beta ~y` performs\n",
    "$\\mathbf{y} \\gets \\alpha \\mathbf{A} \\mathbf{x} + \\beta \\mathbf{y}$ ($\\mathbf{y}$ is destructively modified).\n",
    "Using `gemv`, we write a prediction function of single-layer neural network:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val sigmoid : float -> float = <fun>\n"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "val predict1 :\n",
       "  ('a, 'b, 'c) Slap.D.mat ->\n",
       "  ('a, 'd) Slap_D.vec -> ('b, 'e) Slap.D.vec -> ('a, 'f) Slap_D.vec = <fun>\n"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let sigmoid a = 1.0 /. (1.0 +. exp (~-. a)) ;;\n",
    "let predict1 w b x = Vec.map sigmoid (gemv ~trans:normal w x ~beta:1.0 ~y:(Vec.copy b)) ;;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Two-layer neural networks\n",
    "\n",
    "Let $\\mathbf{x}$ be an input vector of a layer,\n",
    "$\\mathbf{W}_1, \\mathbf{W}_2$ be weight matrices,\n",
    "and $\\mathbf{b}_1, \\mathbf{b}_2$ be biases.\n",
    "Output vector $\\mathbf{z}$ is given by\n",
    "\n",
    "\\begin{align*}\n",
    "  \\mathbf{y} & = \\sigma\\left( \\mathbf{W}_1 \\mathbf{x} + \\mathbf{b}_1 \\right), \\\\\n",
    "  \\mathbf{z} & = \\sigma\\left( \\mathbf{W}_2 \\mathbf{y} + \\mathbf{b}_2 \\right) \\\\\n",
    "\\end{align*}\n",
    "\n",
    "where $\\mathbf{y}$ is a signal of a hidden layer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val predict :\n",
       "  ('a, 'b, 'c) Slap.D.mat ->\n",
       "  ('d, 'a, 'e) Slap.D.mat ->\n",
       "  ('a, 'f) Slap_D.vec ->\n",
       "  ('d, 'g) Slap_D.vec -> ('b, 'h) Slap.D.vec -> ('d, 'i) Slap_D.vec = <fun>\n"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let predict w1 w2 b1 b2 x = predict1 w2 b2 (predict1 w1 b1 x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Training\n",
    "\n",
    "We minimize the root-mean-square (RMS) error\n",
    "\n",
    "$$E(S) = \\frac{1}{2} \\sum_{(\\mathbf{x}, \\mathbf{t}) \\in S} \\| \\mathbf{t} - \\mathbf{z}(\\mathbf{x}) \\|^2$$\n",
    "\n",
    "on a training set $S$ by backpropagation [Rumelhart, et al. 1986]: we repeatedly update the parameters as follows.\n",
    "\n",
    "1. Obtain $(\\mathbf{x}, \\mathbf{t}) \\in S$ with replacement.\n",
    "2. $\\mathbf{y} = \\sigma\\left( \\mathbf{W}_1 \\mathbf{x} + \\mathbf{b}_1 \\right)$\n",
    "3. $\\mathbf{z} = \\sigma\\left( \\mathbf{W}_2 \\mathbf{y} + \\mathbf{b}_2 \\right)$\n",
    "4. $\\mathbf{\\delta}_2 = (\\mathbf{z} - \\mathbf{t}) \\otimes \\mathbf{z} \\otimes (\\mathbf{1} - \\mathbf{z})$\n",
    "5. $\\mathbf{\\delta}_1 = \\mathbf{W}_2^\\top \\mathbf{\\delta}_2 \\otimes \\mathbf{y} \\otimes (\\mathbf{1} - \\mathbf{y})$\n",
    "6. $\\mathbf{W}_1 \\gets \\mathbf{W}_1 - \\eta \\mathbf{\\delta}_1 \\mathbf{x}^\\top$\n",
    "7. $\\mathbf{W}_2 \\gets \\mathbf{W}_2 - \\eta \\mathbf{\\delta}_2 \\mathbf{y}^\\top$\n",
    "8. $\\mathbf{b}_1 \\gets \\mathbf{b}_1 - \\eta \\mathbf{\\delta}_1$\n",
    "9. $\\mathbf{b}_2 \\gets \\mathbf{b}_2 - \\eta \\mathbf{\\delta}_2$\n",
    "\n",
    "where $\\eta > 0$ is a learning rate, and $\\otimes$ is element-wise multiplication on vectors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val train :\n",
       "  eta:float ->\n",
       "  ('a, 'b, 'c) Slap.D.mat ->\n",
       "  ('d, 'a, 'e) Slap.D.mat ->\n",
       "  ('a, 'f) Slap_D.vec ->\n",
       "  ('d, 'g) Slap_D.vec -> ('b, 'h) Slap.D.vec -> ('d, 'i) Slap_D.vec -> unit =\n",
       "  <fun>\n"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let train ~eta w1 w2 b1 b2 x t =\n",
    "  let dsigmoid v =\n",
    "    let ones = Vec.make1 (Vec.dim v) in\n",
    "    Vec.mul v (Vec.sub ones v)\n",
    "  in\n",
    "  let y = Vec.map sigmoid (gemv ~trans:normal w1 x ~beta:1.0 ~y:(Vec.copy b1)) in\n",
    "  let z = Vec.map sigmoid (gemv ~trans:normal w2 y ~beta:1.0 ~y:(Vec.copy b2)) in\n",
    "  let delta2 = Vec.mul (Vec.sub z t) (dsigmoid z) in\n",
    "  let delta1 = Vec.mul (gemv ~trans:trans w2 delta2) (dsigmoid y) in\n",
    "  ignore (ger ~alpha:(~-. eta) delta2 y w2) ;\n",
    "  ignore (ger ~alpha:(~-. eta) delta1 x w1) ;\n",
    "  axpy ~alpha:(~-. eta) b2 delta2 ;\n",
    "  axpy ~alpha:(~-. eta) b1 delta1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Training a random dataset\n",
    "\n",
    "We use a dataset generated from a bivariate Gaussian distribution\n",
    "by [random_dataset_generation.ipynb](random_dataset_generation.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val dataset :\n",
       "  ((Slap.Size.two, 'a) Slap.D.vec * (Slap.Size.one, 'b) Slap.D.vec)\n",
       "  Core.List.t =\n",
       "  [(       R1      R2\n",
       "    -0.740868 1.08918, R1\n",
       "                        0);\n",
       "   (      R1        R2\n",
       "    -1.37741 -0.485505, R1\n",
       "                         0);\n",
       "   (       R1       R2\n",
       "    -0.727394 0.973983, R1\n",
       "                         0);\n",
       "   (       R1      R2\n",
       "    -0.784392 2.71865, R1\n",
       "                        0);\n",
       "   (         R1      R2\n",
       "    -0.00955137 3.33509, R1\n",
       "                          0);\n",
       "   (     R1       R2\n",
       "    0.34222 0.660543, R1\n",
       "                       0);\n",
       "   (       R1        R2\n",
       "    -0.722362 -0.405427, R1\n",
       "                          0);\n",
       "   (      R1       R2\n",
       "    0.146598 -0.49812, R1\n",
       "                        0);\n",
       "   (       R1       R2\n",
       "    -0.811723 -1.19668, R1\n",
       "                         0);\n",
       "   (...); ...]\n"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let dataset =\n",
    "  let float = float_of_string in\n",
    "  In_channel.read_lines \"datasets/bivariate_gaussian_2d.csv\"\n",
    "  |> List.filter_map ~f:(fun line ->\n",
    "      match String.split ~on:',' line with\n",
    "      | [x; y; t] -> Some ([%vec [float x; float y]], [%vec [float t]])\n",
    "      | _ -> None) ;;"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "val w1 : (Slap.Size.five, Slap.Size.two, 'a) Slap_D.mat =\n",
       "             C1        C2\n",
       "  R1     2.3058 -0.107382\n",
       "  R2   -3.27364  0.202057\n",
       "  R3   -3.79329   3.95299\n",
       "  R4    3.15454 -0.161077\n",
       "  R5 -0.0032972 0.0913999\n",
       "val b1 : (Slap.Size.five, 'a) Slap_D.vec =\n",
       "        R1        R2       R3         R4        R5\n",
       "  -0.32692 -0.568707 0.773924 -0.0921993 -0.730415\n",
       "val w2 : (Slap.Size.one, Slap.Size.five, 'a) Slap_D.mat =\n",
       "          C1       C2       C3      C4       C5\n",
       "  R1 2.51502 -2.84573 -2.06153 3.77482 0.566059\n",
       "val b2 : (Slap.Size.one, 'a) Slap_D.vec =      R1\n",
       "                                          0.18147\n"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let (w1, b1, w2, b2) =\n",
    "  let w1 = Mat.random hidden_dim input_dim in\n",
    "  let b1 = Vec.random hidden_dim in\n",
    "  let w2 = Mat.random output_dim hidden_dim in\n",
    "  let b2 = Vec.random output_dim in\n",
    "  for i = 0 to 100 do\n",
    "    List.iter ~f:(fun (x, t) -> train ~eta:1.0 w1 w2 b1 b2 x t) dataset\n",
    "  done ;\n",
    "  (w1, b1, w2, b2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Decision boundary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfQAAAEsCAIAAAC62dafAAAABmJLR0QA/wD/AP+gvaeTAAAYVElEQVR4nO2d0WtU16KH11yszTlJ0YRUfDEtKcVqsa+ObcUglAs9BSHWGhDT6lNIkfjmQwo++Sa3KW1CkHqK4YYSY1MOGFulOqGinf4BFwqFFPLQvIhjHWcwmGHfh5FpdoyzdyZr9vrNmu97SubT5GfndC8z+zvTVBAEBgAA/OK/XA8AAAD7cHEHAPAQLu4AAB7CxR0AwEO4uAMAeAgXdwAAD+HiDgDgIVzcAQA8hIs7AICHcHEHAPAQLu4AAB7CxR0AwEO4uAMAeAgXdwAAD+HiDgDgIVzcAQA8hIs7AICHcHEHAPAQLu4AAB7CxR0AwEO4uAMAeAgXdwAAD+HiDgDgIVzcAQA8hIs7AICHCF3c33zTzM2FHpmbM2++icU2nRWdhV3HUyhAIMPt28G2bcHNm08//fnn0KdYbPNY0VnYdTyF7hG6uAcr/nGt+Q8Ki20eKzoLG9+6RuviHgRBJhO8/HLw8stBJtMstr39SXv7E7VVWFu25udX94+EXWGftLc/aW9/zlPokk2uXxZamwcPHly4cPnHH+dXPb6w0F0oHDPGXLgw6Y3N548YYy5cmJZahbVla35+uxcWjhUKxpjJCxfmf/wRq2kP//VXS0uL4pXU9ekSovLDzYkTF5/9KadiV73Y1ei2t/erI0e+UFuFtWVre36l/0jYFfbTXbsunjjByzLVePXVYOvWp/98zp49e/t2sHVr8Mora9ggCHyyQ0NDQ0NDaquwtmwNz6/70djYdu/evWfPnn3mKXSPUAq5iiDAYrFOvzHWlnWC69MlROVHnxMn/q1z07ve9vjxbz7++Gu1VVhbtrbnV/qPhF1hT7/11uVTp3hZJpryrel//vNRHW5ri9qOjuWOjmW1VVhbtubnV/ePhF1h8//4x/OfQpfoviwDAAC14/p0CdGctczk5OL09KLaKqwtW9vzK/1Hwq6wX/X2PucpdIzQxb1pa5k//vjjjz/+UFuFtWVreH7dj8bGtkNDQ2tYAXRfltG86Y3FJm9FZ2HjWye4Pl1CNGctMzm5ePnyGj+2K2/Gxre1Pb/SfyTsCvtVb+9znkLHaF3cg4BaRmUV1pallvHbUssAAECCuD5dQlDL6KzC2rLUMn5baploqGWkVmFtWWoZvy21zLrRvOmNxSZvRWdh41snuD5dQlDL6KzC2rLUMn5bapm4ZKhlNFZhbVlqGb8ttQwAACRIwofJDz/8cPDgwZaWlh07dpw+ffrRo0crLbWMziqsLUst47ellnnKe++9d+3atXw+/+eff/b19Z08ebKiqGWkVmFtWWoZv61sLZP0f9b1xo0b5Q/a2tpGR0dff/31ixcvrvkrNW96Y7HJW9FZ2PjWCQ4Plt9++627u3vlI9QyOquwtiy1jN9WtpZJBe4OnL6+vj179gwPD698cG7O9PYuF4vFw4enurrmjTGtra39/f1lm822DA52GmPGxu6l049XfcEGtQMDHcaY8fH7UquwtmzNz29LNts5OGiMuTc29jidxurYiYmJQqFgjOleWDgyPd3W1rZpZsb09Bglkn5ZpsLIyEgulztz5syaNgiCfD6fy+WMMalUanl5ufx4qVQqn0alUqnyYIUGtZUPpFZhbdnKB1KrsBuxxWKxfHXK5/MO/34cQf1+KKjC+fPn9+/fXygUVj1OLaOzCmvLUsv4ball/ubSpUvpdPrhw4erHqeWkVqFtWWpZfy21DJP+emnn0ZHR69fv/7SSy9V/5WaN72x2OSt6CxsfOuEhA+TLVu2rBqQy+UqllpGZxXWlqWW8dvK1jK8t4x7y3vL+G15bxm/Le8tAwAACeL6dAlBLaOzCmvLUsv4balloqGWkVqFtWWpZfy2srWM7ssymje9sdjkregsbHzrBNenSwhqGZ1VWFuWWsZvSy0Tlwy1jMYqrC1LLeO3pZYBAIAEcX26hKCW0VmFtWWpZfy21DLRUMtIrcLastQyfltqmXWjedMbi03eis7CxrdOcH26hKCW0VmFtWWpZfy21DJxyVDLaKzC2rLUMn5bahkAAEgQ16dLCGoZnVVYW5Zaxm9LLRMNtYzUKqwtu+bz+/Orx7dtXVr5e7dtXbr5ygmV0djYllpm3Wje9MZirdgBMz5ljh40t8qfvhvcnjJHT5kvHc/C1sk6wfXpEoJaRmcV1pZ93vMrPRob21LLxCVDLaOxCmvLVnl+dUdjY1tqGQAASBDXp0sIahmdVVhb9nnPr/RobGxLLRMNtYzUKqwtu+bz634W1pKlllk3mje9sdh6WNFZWFvWCa5PlxDUMjqrsLYstYzfllomLhlqGY1VWFuWWsZvSy0DAAAJ4vp0CUEto7MKa8tSy/htqWWioZaRWoW1Zall/LbUMutG86Y3FlsPKzoLa8s6wfXpEoJaRmcV1pallvHbUsvEJUMto7EKa8tSy/htqWUAACBBXJ8uIahldFZhbVlqGb8ttUw01DJSq7C2LLWM35ZaZt1o3vTGYuthRWdhbVknuD5dQlDL6KzC2rLUMn5bapm4ZKhlNFZhbVlqGb8ttQwAACSIkyPl5s2b6XT62e9OLaOzCmvLUsv4ballQhw4cCCTyay6uFPLSK3C2rLUMn5b2Vpmk5MfF+bm5iJ/jeZNbyy2HlZ0FtaWdYLDg+XZ704to7MKa8tSy/htZWuZVODuwEml1vjuc3Omt3e5WCwePjzV1TVvjGltbe3v7y/bbLZlcLDTGDM2di+dfrzq9zaoHRjoMMaMj9+XWoW1Zas8vy3ZbOfgoDHm3tjY43Qa2yh2YmKiUCgYY7oXFo5MT7e1tW2amTE9PUYJNy/LRBIEQT6fz+VyxphUKrW8vFx+vFQqlc+DUqlUebBCg9rKB1KrsLZs5QOpVdiN2GKxWL465fN5h38/jqBuPxNE8+x3p5bRWYW1Zall/LbUMmt9b2oZapkmsNQyflvZWkb3/8SkedMbi62HFZ2FtWWd4ORIed4GahmdVVhbllrGbytby/DeMu4t7y3jt+W9Zfy2vLcMAAAkiOvTJQS1jM4qrC1LLeO3pZaJhlpGahXWlqWW8dtSy6wbzZveWGw9rOgsrC3rBNenSwhqGZ1VWFuWWsZvSy0Tlwy1jMYqrC1LLeO3pZYBAIAEcX26hKCW0VmFtWWpZfy21DLRUMtIrcLastQyfltqmXWjedMbi62HFZ2FtWWd4Pp0CUEto7MKa8tSy/htqWXikqGW0ViFtWWpZfy21DIAAJAgrk+XENQyOquwtiy1jN+WWiYaahmpVdgN2p9fPb5t69LK53fb1qWbr5yQHo2llkkAzZveWGxMO2DGp8zRg+ZW+dN3g9tT5ugp86XjWVgn1gmuT5cQ1DI6q7Abt2W9ODm5ePmy0CysVUstE5cMtYzGKqwVG2Qyyx0dyx0dWrOw9iy1DAAAJIjr0yUEtYzOKuzGbVkvTk4uTk8LzcJatdQy0VDLSK3CbtBWdPn5VZmFpZZxjuZNbyy2Nis6C5uMdYLr0yUEtYzOKuzGLbVMM1hqmbhkqGU0VmGtWGoZ7y21DAAAJIjr0yUEtYzOKuzGLbVMM1hqmWioZaRWYTdoqWWaxFLLrBvNm95YbG1WdBY2GesE16dLCGoZnVXYjVtqmWaw1DJxyVDLaKzCWrHUMt5bahkAAEgQ16dLCGoZnVXYjVtqmWaw1DLRUMtIrcJG2sp/a6li1/xvLVHL+G2pZdaN5k1vLLYC/60lbFzrBNenSwhqGZ1V2Dg2zm+mlvHbUsvEJUMto7EKG9NG/mZqGb8ttQwAACRI8ufJlStXdu3atXnz5l27ds3MzKxU1DI6q7BxbJzfTC3jt6WWecrdu3c7OztnZ2cfPHgwOzvb2dmZzWbLilpGahU20sb8zdQyfltqmaeMjIwMDw+///77W7Zsef/994eHh0dGRp7zI0W1r4PFqlnRWVgF64SED5MdO3b8/vvvlU9///33rq6uyqfUMjqrsHEstQxWtpZJBckeOJs3b87lcq2treVPC4VCR0fH0tJS5RfMzZne3uVisXj48FRX17wxprW1tb+/v2yz2ZbBwU5jzNjYvXT68aov3qB2YKDDGDM+fl9qFTambclmOwcHjTH3xsYep9PP2o6BAWPM/fHxNW3134vVtBMTE4VCwRjTvbBwZHq6ra1t08yM6ekxSmxyPWBtgiDI5/O5XM4Yk0qllpeXy4+XSqXyaVQqlSoPVmhQW/lAahXWlq18ILUKuxFbLBbLV6d8Pp/w34/XQf1+KFiT6i/LUMvorMLGsdQyWGqZp3z00Ueff/555dPPP/+8r6+v/DG1jNQqbBD73WOqf2lqGb+tbC2T9MsyQ0NDhw4d2rlz5zvvvHPnzp1z585dvXr1OT9SVPs6WGwCtvzuMT3mlDEHzdN3j/mfT83//p/bWdiGs05I/jyZnp7euXPnCy+88MYbb3z33XcrFbWMziqsrS9NLeO3la1leG8Z95b3lhG3G/zSvLeM35b3lgEAgARxfbqEoJbRWYW19aWpZfy21DLRUMtIrcLa+tLUMn5b2VpG92UZzZve2Ga2orOw+tYJrk+XENQyOquwtr40tYzfllomLhlqGY1VWFtfmlrGb0stAwAACeL6dAlBLaOzCmvrS1PL+G2pZaKhlpFa1STWyrvHVLfUMn5bapl1o3nTG+uZLb97zEFzq/zpu8HtKXP0lPnS8SysZ9YJrk+XENQyOquax9b7G1PL+G2pZeKSoZbRWNVUtq7fmFrGb0stAwAACeL6dAlBLaOzqnlsvb8xtYzfllomGmoZqVVNYhP4xtQyfltqmXWjedMb67cVnYVtdOsE16dLCGoZnVXNY6llsBux1DJxyVDLaKxqKkstg6WWAQCARsD16RKCWkZnVfNYahnsRiy1TDTUMlKrvLEJvHtMdUst47elllk3mje9sQ1nefcYrHvrBNenSwhqGZ1VPlm3s6hl/LbUMnHJUMtorPLMOpxFLeO3pZYBAIAEcX26hKCW0Vnlk3U7i1rGb0stEw21jNQqb6zzWdQyfltqmXWjedMb2+hWdBbWb+sE16dLCGoZnVU+WWoZbP0stUxcMtQyGqs8s9QyWGoZAABofFyfLiGoZXRW+WSpZbD1s9Qy0VDLSK1qIOv83WOqW2oZvy21zLrRvOmNFbS8ewxW3TrB9ekSglpGZ1VjWdFZ1DJNYKll4pKhltFY1XBWdBa1jO+WWgYAABLEyZFy8+bNdDr97HenltFZ1VhWdBa1TBNYapkQBw4cyGQyqy7u1DJSqxrIis6ilmkOK1vLbHLy48Lc3Fzkr9G86Y3Vt6KzsM1sneDwYHn2u1PL6KxqLCs6i1qmCaxsLZMK3B04qdQa331uzvT2LheLxcOHp7q65o0xra2t/f39ZZvNtgwOdhpjxsbupdOPV/3eBrUDAx3GmPHx+1KrGs62ZLOdg4PGmHtjY4/TaR3bMTBgjLk/Pi61CrsROzExUSgUjDHdCwtHpqfb2to2zcyYnh6jhJuXZSIJgiCfz+dyOWNMKpVaXl4uP14qlcrnQalUqjxYoUFt5QOpVVhbtvKB1CrsRmyxWCxfnfL5vMO/H0dQt58J/uZ53+7Z704to7OqsazoLGqZJrDUMmt9b2oZapl4VvzdY6pbahm/rWwto/t/YtK86Y11Ynn3GGxjWyc4OVKet4FaRmeVmhWdRS3T9Fa2luG9Zdxb3lsmphWdFWV5bxm/Le8tAwAACeL6dAlBLaOzSs2KzqKWaXpLLRMNtYzUKikrOotaBkstUwOaN72xClZ0FhYrhevTJQS1jM4qNSs6i1qm6S21TFwy1DIaqwSt6KwoSy3jt6WWAQCABHF9uoSgltFZpWZFZ1HLNL2llomGWkZqlZQVnUUtg6WWqQHNm95YBSs6C4uVwvXpEoJaRmeVmhWdRS3T9JZaJi4ZahmNVYJWdFaUpZbx21LLAABAgrg+XUJQy+isUrOis6hlmt5Sy0RDLSO1SsqKzqKWwVLL1IDmTW+sghWdhcVK4fp0CUEto7NKzYrOopZpekstE5cMtYzGKkErOivKUsv4ballAAAgQVyfLiGoZXRWqVnRWdQyTW+pZaKhlpFaJWVFZ1HLYKllakDzpjdWwYrOwmKlcH26hKCW0VmlZkVnUcs0vaWWiUuGWkZjlaAVnRVlqWX8ttQyAACQIK5PlxDUMjqr1KzoLGqZprfUMtFQy0itStj+/OrxbVuXVtptW5duvnJCejS1DJZapgY0b3pj62QHzPiUOXrQ3Cp/+m5we8ocPWW+dDwLi7VineD6dAlBLaOzKnkrOotaBlvVUsvEJUMto7GKHsaWpZbx21LLAABAgrg+XUJQy+isSt6KzqKWwVa11DLRUMtIrUrYis6ilsFGWWqZdaN50xubjBWdhcXWZp3g+nQJQS2jsyp5KzqLWgZb1VLLxCVDLaOxilrGlqWW8dtSywAAQIIkfJj88MMPBw8ebGlp2bFjx+nTpx89erTSUsvorEreis6ilsFWtdQyT3nvvfeuXbuWz+f//PPPvr6+kydPVhS1jNSqhK3oLGoZbJSVrWU2JfyDwo0bN8oftLW1jY6Ovv766xcvXlzzV2re9MYmY0VnYbG1WSc4PFh+++237u7ulY9Qy+isSt6KzqKWwVa1srVMKnB34PT19e3Zs2d4eHjlg3Nzprd3uVgsHj481dU1b4xpbW3t7+8v22y2ZXCw0xgzNnYvnX686gs2qB0Y6DDGjI/fl1rlxLZks52Dg8aYe2Njj9NpP2zHwIAx5v74uNQq7EbsxMREoVAwxnQvLByZnm5ra9s0M2N6eowSSb8sU2FkZCSXy505c2ZNGwRBPp/P5XLGmFQqtby8XH68VCqVT6NSqVR5sEKD2soHUquwtmzlA6lV2I3YYrFYvjrl83mHfz+OoH4/FJRZ83udP39+//79hUJh1S+mltFZlbwVnUUtg61qqWX+5tKlS+l0+uHDh6sep5aRWpWwFZ1FLYONstQyT/npp59GR0evX7/+0ksvVf+Vmje9sclY0VlYbG3WCQkfJlu2bFk1IJfLVSy1jM6q5K3oLGoZbFUrW8vw3jLuLe8t4/4b19Py3jJ+W95bBgAAEsT16RKCWkZnVfJWdBa1DLaqpZaJhlpGalXCVnQWtQw2ysrWMrovy2je9MYmY0VnYbG1WSe4Pl1CUMvorEreis6ilsFWtdQycclQy2isopaxZall/LbUMgAAkCCuT5cQ1DI6q5K3orOoZbBVLbVMNNQyUqsStqKzqGWwUZZaZt1o3vTGJmNFZ2GxtVknuD5dQlDL6KxK3orOopbBVrXUMnHJUMtorKKWsWWpZfy21DIAAJAgrk+XEM1Zyxw//s0nn3yttip5Kzprw/ab48e//uQTtVVYW/b0W29NnTrFyzLVaNpaZmhoaGhoSG1VwlZ0lqWaYmhoSG0V1pbdu3fv2bNnqWXWgeZNb2wyVnQWFlubdYLr0yVEc9Yyvb1fffjhF2qrkreis2zUFF98+KHaKqwt++muXf8+cULwZZlUIHbgzM2Zjz4yf/311yefXOvuXlhl5+e7JicPGWOOHfuPN/bSpX8ZYz7+eFZqlRPbNT9/aHLSGPOfY8cWurv9sP+6dMkYM/vxx1KrsLbsf1+8+OKLL7bOzpqeHiOF69NlNfVslkRte/uT9vYnaqtIIW3ZJ+3tT9rb1VZhbdlqz69TtC7uwrfEsXW3orOw2Br/B+0YoYv7qhe1Vv3jwvptRWdhsTX+D9o9Qhf33btX/2STyQS7d2ObworOwmJrswLI3VAFAICNo9u5AwBAzXBxV+HWrVv79u1LpVKuh4BNvvvuu927d7/44ou7d+/+/vvvXc8By0j/a+v6dSF4yoEDBzKZDM+IT9y9e7ezs3N2dvbBgwezs7OdnZ3ZbNb1KLCJ8r+2vOauRSrFM+IPR48e3bdv3+nTp8ufjoyM/Prrr99++63bVWAdzX9teVkGoF788ssvH3zwQeXTDz744O7duw73QFOheOA0M5p/BYDa2Lx5cy6Xa21tLX9aKBQ6OjqWlpbcrgLraP5ry9/cAQA8hIu7G1IrcL0F6sX27dsXFxcrny4uLm7fvt3hHmgquLi7YeVNbddboF7s27fv6tWrlU+vXr369ttvO9wDTcUm1wMAvGVoaOjQoUM7d+5855137ty5c+7cuZXXeoC6ongfoDlZ9foMz4sfXLly5bPPPpufn3/ttdfOnTvX29vrehHYRPlfWy7uAAAewmvuAAAewsUdAMBDuLgDAHgIF3cAAA/h4g4A4CFc3AEAPISLOwCAh3BxBwDwEC7uAAAewsUdAMBDuLgDAHgIF3cAAA/h4g4A4CH/D06q7wUxCBP/AAAAAElFTkSuQmCC"
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "- : unit = ()\n"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "let xmin, xmax = -1.5, 1.0 in\n",
    "let ymin, ymax = -2.5, 2.5 in\n",
    "let n_samples = 40 in\n",
    "let vp = A.init ~text:14. ~w:500. ~h:300. [\"jupyter\"] in\n",
    "A.Axes.box vp ;\n",
    "A.Viewport.xrange vp xmin xmax ;\n",
    "A.Viewport.yrange vp ymin ymax ;\n",
    "for i = 0 to n_samples do\n",
    "  let y = float i /. float n_samples *. (ymax -. ymin) +. ymin in\n",
    "  for j = 0 to n_samples do\n",
    "    let x = float j /. float n_samples *. (xmax -. xmin) +. xmin in\n",
    "    let z = predict w1 w2 b1 b2 @@ Vec.of_list_dyn input_dim [x; y] in\n",
    "    A.set_color vp (if Vec.get_dyn z 1 > 0.5 then A.Color.red else A.Color.blue) ;\n",
    "    A.List.xy_pairs vp ~style:(`Markers \"x\") [x, y]\n",
    "  done\n",
    "done ;\n",
    "A.close vp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "OCaml 4.04.2",
   "language": "OCaml",
   "name": "ocaml-jupyter-4.04.2"
  },
  "language_info": {
   "codemirror_mode": "text/x-ocaml",
   "file_extension": ".ml",
   "mimetype": "text/x-ocaml",
   "name": "OCaml",
   "nbconverter_exporter": null,
   "pygments_lexer": "OCaml",
   "version": "4.04.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}