sparse_latent_semantic_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import nltk\n",
    "from nltk.corpus import stopwords\n",
    "import itertools\n",
    "from collections import Counter\n",
    "from sklearn.decomposition import PCA\n",
    "from sklearn.metrics.pairwise import cosine_similarity\n",
    "from tqdm import tqdm\n",
    "\n",
    "nltk.download('brown')\n",
    "nltk.download('stopwords')\n",
    "stopwords = stopwords.words('english')\n",
    "\n",
    "np.set_printoptions(suppress=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Latent Semantic Analysis (LSA)\n",
    "\n",
    "Word representations can be either statistics-based or learning based (Word2Vec, ELMo, BERT embeddings and more..). In a statistics based approach, we use the global count statistics on a large corpus to serve as a starting point to model semantic and syntactic relationships. Since the corpus can be very large comprising of large number documents and rich vocabulary, encapsulating this information in a structured format as matrix and further processing it can be memory inefficient and computationally costly.\n",
    "\n",
    "The idea of LSA is to encode word relationships in lower dimensions. LSA essentially performs **Singular Value Decomposition (SVD)** on the term-document matrix reducing the number of dimensions to a finite number (say number of topics) - which is where it is used for topic modelling. This enables us to represent documents from the original word space (every document (row in term-document matrix) can be uniquely identified by the terms it contains) to a latent topic space (document represented by a continuous/weighted distribution of topics). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# info\n",
    "\n",
    "# number of training epochs for Sparse LSA\n",
    "n_epochs = 4000\n",
    "\n",
    "# number of documents to consider\n",
    "n_docs = 10\n",
    "\n",
    "# number of topics\n",
    "n_topics = 5\n",
    "\n",
    "seed = 42"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Processing sentences..\n",
      "\n",
      "Number of Documents: 10\n",
      "Number of Tokens: 4447\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(seed)\n",
    "brown = nltk.corpus.brown\n",
    "docs = np.random.choice(brown.fileids(), n_docs, replace=False)\n",
    "\n",
    "print('Processing sentences..\\n')\n",
    "processed_docs = {}\n",
    "for doc in docs:\n",
    "    processed_sents = []\n",
    "    sents = brown.sents(doc)\n",
    "    for sent in sents:\n",
    "        processed_sents.append([word.lower() for word in sent if word.isalnum() and word not in stopwords])\n",
    "    processed_docs[doc] = processed_sents\n",
    "\n",
    "processed_sents = list(itertools.chain(*list(processed_docs.values())))\n",
    "final_tokens = list(set(itertools.chain(*processed_sents)))\n",
    "n_tokens = len(final_tokens)\n",
    "\n",
    "token2int = dict(zip(final_tokens, range(n_tokens)))\n",
    "int2token = {v:k for k,v in token2int.items()}\n",
    "\n",
    "print('Number of Documents:', n_docs) \n",
    "print('Number of Tokens:', n_tokens)    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_frequency_hash_map(tokens):\n",
    "    token_frequency_dict = {}\n",
    "    all_tokens = list(itertools.chain(*tokens))\n",
    "    for token in all_tokens:\n",
    "        if token_frequency_dict.get(token):\n",
    "            token_frequency_dict[token] +=1\n",
    "        else:\n",
    "            token_frequency_dict[token] = 1\n",
    "    return token_frequency_dict"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def build_term_document_matrix(processed_docs):\n",
    "    term_document_matrix = np.zeros((n_docs, n_tokens), dtype='int')\n",
    "    for doc_idx, doc in enumerate(processed_docs):\n",
    "        tokens = processed_docs[doc]\n",
    "        token_frequency_dict = get_frequency_hash_map(tokens)\n",
    "        token_idxs = list(map(lambda x: token2int[x], token_frequency_dict.keys()))\n",
    "        counts = list(token_frequency_dict.values())        \n",
    "        term_document_matrix[doc_idx, token_idxs] = counts\n",
    "    assert term_document_matrix.shape == (n_docs, n_tokens)        \n",
    "    return term_document_matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_singular_value_decomposition(term_document_matrix, n_topics):\n",
    "    U, S, VT = np.linalg.svd(term_document_matrix)\n",
    "    \n",
    "    \"\"\"\n",
    "    Necessary checks for Singular Value Decomposition:\n",
    "\n",
    "    1. U and V are singular matrices\n",
    "    2. Number of eigenvalues is rank of matrix \n",
    "    \"\"\"\n",
    "\n",
    "    assert np.allclose(U.T @ U, np.identity(len(U)))\n",
    "    assert np.allclose(VT @ VT.T, np.identity(len(VT)))\n",
    "    assert len(S) == np.linalg.matrix_rank(term_document_matrix)\n",
    "    \n",
    "    # low-rank approximation\n",
    "    approximated_mat = U[:,:n_topics] @ np.diag(S[:n_topics]) @ VT[:n_topics]\n",
    "    assert approximated_mat.shape == term_document_matrix.shape\n",
    "\n",
    "    # projection matrix\n",
    "    projection_mat = np.linalg.inv(np.diag(S[:n_topics])) @ VT[:n_topics]\n",
    "    assert projection_mat.shape == (n_topics, n_tokens)\n",
    "    return approximated_mat, projection_mat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Building Term Document Matrix..\n",
      "Getting document representations in topic space (instead of original word space)..\n"
     ]
    }
   ],
   "source": [
    "print('Building Term Document Matrix..')\n",
    "term_document_matrix = build_term_document_matrix(processed_docs)\n",
    "\n",
    "print('Getting document representations in topic space (instead of original word space)..')\n",
    "approximated_mat, projection_mat = get_singular_value_decomposition(term_document_matrix, n_topics)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's look at similarities between document representations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "csim = cosine_similarity(projection_mat)\n",
    "assert np.allclose(csim, np.identity(n_topics))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This says that the document representations in latent topic space i.e. **topics embeddings** are orthogonal to each other. This means that the original documents selected comprised of topics different from each other. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sparse LSA\n",
    "\n",
    "From the [paper](https://www.cs.cmu.edu/~jgc/publication/PublicationPDF/Sparse_Latent_Semantic_Analysis.pdf), we are interested in solving the optimization problem:\n",
    "\n",
    "\\begin{aligned}\n",
    "J = \\min _{\\mathbf{U}, \\mathbf{A}} & \\frac{1}{2}\\|\\mathbf{X}-\\mathbf{U} \\mathbf{A}\\|_{F}^{2}+\\lambda\\|\\mathbf{A}\\|_{1} \\\\\n",
    "\\text { subject to: } & \\mathbf{U}^{T} \\mathbf{U}=\\mathbf{I}\n",
    "\\end{aligned}\n",
    "\n",
    "Differentiating the loss function $J$ with respect to $U$ and $A$:\n",
    "\n",
    "\\begin{aligned}\n",
    "&\\frac{\\partial J}{\\partial A}=A-U^{\\top} X+\\lambda \\operatorname{sgn}(A) \\\\\n",
    "\n",
    "&\\frac{\\partial J}{\\partial U}=2 \\lambda U-X A^{\\top}\n",
    "\\end{aligned}\n",
    "\n",
    "Note that the loss is scaler while $U$ and $A$ are matrices and so we are differentiating a scaler with a matrix."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Derivation \n",
    "\n",
    "We have the governing equation:\n",
    "\n",
    "\\begin{aligned}\n",
    "&X=U A+\\epsilon\n",
    "\\end{aligned}\n",
    "\n",
    "The loss $J$ is given by sum of squares of errors:\n",
    "\\begin{aligned}\n",
    "&J=\\frac{\\|\\epsilon\\|^{2}}{2}=\\frac{\\operatorname{tr}\\left(\\epsilon^{\\top} \\epsilon\\right)}{2}=\\frac{\\operatorname{tr}\\left((X-U A)^{\\top}(X-U A)\\right)}{2}\n",
    "\\end{aligned}\n",
    "\n",
    "$$\n",
    "\\frac{\\operatorname{tr}\\left((X-U A)^{\\top}(X-U A)\\right)}{2}=\\frac{\\operatorname{tr}\\left(X^{\\top} X-X^{\\top} U A-A^{\\top} U^{\\top} X+A^{\\top} U^{\\top} U A\\right)}{2}\n",
    "$$\n",
    "\n",
    "Differentiating with respect to $U$ and $A$:\n",
    "\\begin{aligned}\n",
    "&\\frac{1}{2} \\frac{\\partial \\operatorname{tr}\\left(A^{\\top} A\\right)}{\\partial A}=A \\\\\n",
    "&\\frac{\\partial}{\\partial A} \\operatorname{tr}\\left(A^{\\top} U^{\\top} X\\right)=U^{\\top} X \\\\\n",
    "&\\frac{\\partial}{\\partial U}\\left(U^{\\top} U\\right)=2 U \\\\\n",
    "&\\frac{\\partial}{\\partial U} \\operatorname{tr}\\left(A^{\\top} U^{\\top} X\\right)=X A^{\\top}\n",
    "\\end{aligned}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def is_present(token, words):\n",
    "    if token in words:\n",
    "        return 1\n",
    "    else:\n",
    "        return 0\n",
    "\n",
    "def get_present_label(doc):\n",
    "    words = brown.words(doc)\n",
    "    words = [word.lower() for word in words if word.isalnum() and word not in stopwords] \n",
    "    present_label = list(map(lambda x: is_present(x, words), final_tokens))\n",
    "    return present_label\n",
    "    \n",
    "def get_term_freq(doc):\n",
    "    words_doc = list(itertools.chain(*processed_docs[doc]))\n",
    "    vocab_doc = set(words_doc)    \n",
    "\n",
    "    term_freq_dict = dict(Counter(words_doc))\n",
    "    term_freq_dict = {token2int[k]: v for k,v in term_freq_dict.items()}\n",
    "    counts = np.zeros(n_tokens)\n",
    "\n",
    "    term_freq_mat = np.array(sorted(list(zip(term_freq_dict.keys(), term_freq_dict.values()))))\n",
    "    idxs = term_freq_mat[:, 0]\n",
    "    values = term_freq_mat[:, 1]\n",
    "    counts[idxs] = values\n",
    "    return counts / len(vocab_doc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# term frequency matrix\n",
    "term_freq_mat = np.array(list(map(lambda x: get_term_freq(x), docs)))\n",
    "\n",
    "# document frequency matrix\n",
    "all_docs = brown.fileids()\n",
    "doc_frequency_mat = np.array(list(map(lambda x: get_present_label(x), all_docs)))\n",
    "idf = np.log(len(all_docs) / np.sum(doc_frequency_mat, axis=0))\n",
    "\n",
    "# tf-idf matrix\n",
    "tfidf_mat = term_freq_mat * idf\n",
    "assert tfidf_mat.shape == (n_docs, n_tokens)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def loss_fn(X,U,A,lamda):\n",
    "    return np.sum((X-U@A).T @ (X-U@A)) / 2 + lamda * np.sum(abs(A))\n",
    "\n",
    "def gradient(X,U,A,reg_coef,lamda):\n",
    "    dA = A - U.T @ X + lamda * np.sign(A)\n",
    "    dU = 2 * reg_coef * U - X @ A.T\n",
    "    return dA, dU\n",
    "\n",
    "def sparse_LSA_optimization(X, alpha, lamda, reg_coef, eps):\n",
    "    U = np.random.random((n_docs, n_topics))\n",
    "    A = np.random.random((n_topics, n_tokens))\n",
    "    dA = np.random.random(A.shape)\n",
    "    dU = np.random.random(U.shape)\n",
    "\n",
    "    costs = []\n",
    "    n_iter = 0\n",
    "    cost = 1\n",
    "    convergence = 1\n",
    "    while cost > eps:\n",
    "        dA, dU = gradient(X,U,A,reg_coef,lamda)\n",
    "        A-= alpha * dA\n",
    "        U-= alpha * dU\n",
    "        cost = loss_fn(X,U,A,lamda)\n",
    "\n",
    "        if n_iter%100==0:\n",
    "            print(f'Cost at {n_iter} iterations:', cost.round(3))\n",
    "        \n",
    "        costs.append(cost)\n",
    "        n_iter += 1\n",
    "\n",
    "        if n_iter>=n_epochs:\n",
    "            convergence = 0\n",
    "            break        \n",
    "\n",
    "    if convergence:\n",
    "        print(f'Converged in {len(costs)} epochs..')\n",
    "    else:\n",
    "        print(f'Training complete with {n_epochs} epochs..')\n",
    "    return A, U, costs        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Cost at 0 iterations: 201152295.776\n",
      "Cost at 100 iterations: 195931038.104\n",
      "Cost at 200 iterations: 160427788.643\n",
      "Cost at 300 iterations: 115621216.588\n",
      "Cost at 400 iterations: 74842211.731\n",
      "Cost at 500 iterations: 43710355.941\n",
      "Cost at 600 iterations: 22933028.313\n",
      "Cost at 700 iterations: 10637132.098\n",
      "Cost at 800 iterations: 4220157.348\n",
      "Cost at 900 iterations: 1343544.548\n",
      "Cost at 1000 iterations: 302405.153\n",
      "Cost at 1100 iterations: 34326.668\n",
      "Cost at 1200 iterations: 325.874\n",
      "Cost at 1300 iterations: 64.827\n",
      "Cost at 1400 iterations: 68.246\n",
      "Cost at 1500 iterations: 70.163\n",
      "Cost at 1600 iterations: 71.173\n",
      "Cost at 1700 iterations: 71.452\n",
      "Cost at 1800 iterations: 71.128\n",
      "Cost at 1900 iterations: 71.143\n",
      "Cost at 2000 iterations: 71.645\n",
      "Cost at 2100 iterations: 71.46\n",
      "Cost at 2200 iterations: 71.849\n",
      "Cost at 2300 iterations: 72.125\n",
      "Cost at 2400 iterations: 72.011\n",
      "Cost at 2500 iterations: 71.733\n",
      "Cost at 2600 iterations: 71.637\n",
      "Cost at 2700 iterations: 71.812\n",
      "Cost at 2800 iterations: 71.737\n",
      "Cost at 2900 iterations: 71.701\n",
      "Cost at 3000 iterations: 71.67\n",
      "Cost at 3100 iterations: 71.612\n",
      "Cost at 3200 iterations: 71.571\n",
      "Cost at 3300 iterations: 71.695\n",
      "Cost at 3400 iterations: 71.768\n",
      "Cost at 3500 iterations: 71.752\n",
      "Cost at 3600 iterations: 71.758\n",
      "Cost at 3700 iterations: 71.773\n",
      "Cost at 3800 iterations: 71.796\n",
      "Cost at 3900 iterations: 71.804\n",
      "Training complete with 4000 epochs..\n"
     ]
    }
   ],
   "source": [
    "X = tfidf_mat\n",
    "alpha = 0.001\n",
    "lamda = 0.4\n",
    "reg_coef = 0.5\n",
    "eps = 0.01    \n",
    "\n",
    "A, U, costs = sparse_LSA_optimization(X, alpha, lamda, reg_coef, eps)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEWCAYAAAB1xKBvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xl8XeV95/HPV6sXyfImW963GK/FBhSzpcSUhjG0KUmHTkwW0jQZlwzklbTTtGTSSdK0M9OmadqkISEOYYAskEwDCcmYpVlYEgJYBgM2xmCMjYVXjPG+Sfr1j3uU3DiSdWXr6Nx79X2/Xvele895zr0/HUv+6pznnOdRRGBmZtaTiqwLMDOz0uDAMDOzgjgwzMysIA4MMzMriAPDzMwK4sAwM7OCODDMzKwgDgyzjEh6QNIHsq7DrFAODDMzK4gDw6wXJE2SdKekXZJ2S/qipApJfy1ps6Sdkm6T1JC0HyTpG0nb1yWtlDRW0v8Cfhv4oqQDyftI0j8n77FX0tOS5mf7HZv9igPDrECSKoEfApuBqcAE4A7gj5PHxcB0oA74YrLZe4EGYBIwCrgGOBwRHwceBq6LiLqIuA64FLgIOAMYDrwD2J3+d2ZWmLILDEk3J3+hrSmg7WRJP5X0ZPLX3OX9UaOVrEXAeOCjEXEwIo5ExM+AdwGfi4iNEXEA+BiwVFIVcJxcULwhItojYlVE7Ovm/Y8D9cBsQBGxLiK2pf5dmRWo7AIDuAVYUmDbvwa+ExFnAUuBL6VVlJWFScDmiGg7Yfl4ckcdnTYDVcBY4OvAfcAdkrZK+oyk6q7ePCJ+Qu7I5AZgh6Tlkob19TdhdqrKLjAi4iHgtfxlkmZIulfSKkkPS5rd2Rzo/IVsALb2Y6lWerYAk5Mjh3xbgSl5rycDbcCOiDgeEX8TEXOBC4DfB65O2v3GUNER8YWIOAeYR+7U1Ef7+HswO2VlFxjdWA58KPlF/At+dSTxKeDdklqBFcCHsinPSsTjwDbg7yUNTTq0LwRuB/5M0jRJdcD/Br4dEW2SLpb0W0n/xz5yp53ak/fbQa7PAwBJb5R0bnIEchA4ktfWLHNlHxjJL/AFwP+TtBr4CjAuWX0VcEtETAQuB74uqez3iZ2aiGgH3gq8AXgZaCXXMX0zuVNPDwEvkfuPvvOPjybg38iFxTrgQeAbybrPA1dK2iPpC+SOdr8K7CF3Wms38NnUvzGzAqkcJ1CSNBX4YUTMT84Br4+IcV20WwssiYgtyeuNwHkRsbM/6zUzKwVl/9d0ckXKS5L+CCC51n1Bsvpl4JJk+RxgELArk0LNzIpc2R1hSLodWAyMJneO+JPAT4AvkzsVVQ3cERGfljSX3CmAOnIdkH8ZEfdnUbeZWbEru8AwM7N0lP0pKTMz6xsnXk9e0kaPHh1Tp07Nugwzs5KxatWqVyOisZC2ZRUYU6dOpaWlJesyzMxKhqTNPbfK8SkpMzMriAPDzMwK4sAwM7OCODDMzKwgDgwzMyuIA8PMzAriwDAzs4I4MPKs3vI6X3pgA49t9DTKZmYnKqsb907Ht1e+zPV3PkPn0Fp/cekZXPc7M7MtysysiKR2hCFpkqSfSlonaa2kD3fRRpK+IGmDpKclnZ23bomk9cm669OqE2DDzgP8z++t5U1vGM3jH7+Et581gc/e/zwPrPe0GGZmndI8JdUG/PeImAOcB1ybDCee7zJgZvJYRm4IcpLpLG9I1s8Frupi2z7zuX9fT3Wl+Od3LGRM/SD+zx/+FjMah/LpHzzL8faOtD7WzKykpBYYEbEtIp5Inu8nNz3lhBOaXQHcFjmPAsMljQMWARsiYmNEHAPuSNr2ub2Hj9OyaQ/vf9M0RtfVAjCoupLrL5vDxlcP8v3VW9P4WDOzktMvnd7JlKlnAY+dsGoCsCXvdWuyrLvlXb33Mkktklp27er9ZHkNg6t58KMX86dvnvFry393zhjeMKaOWx/ZhOcMMTPrh8CQVAd8F/hIMl3qr63uYpM4yfLfXBixPCKaI6K5sbGgEXp/w+CaSobW/nr/vyTee/4UnnllL0+17j2l9zUzKyepBoakanJh8c2IuLOLJq3ApLzXE4GtJ1ner9521gRqqir43pOv9PdHm5kVnTSvkhLwNWBdRHyum2Z3A1cnV0udB+yNiG3ASmCmpGmSaoClSdt+VT+omotnNbLimW20d/i0lJkNbGkeYVwIvAf4HUmrk8flkq6RdE3SZgWwEdgAfBX4bwAR0QZcB9xHrrP8OxGxNsVau/X7Z45n5/6jrNz0WhYfb2ZWNFK7cS8ifkbXfRH5bQK4tpt1K8gFSqYumTOGmsoKfrxuB+dNH5V1OWZmmfHQID0YUlPFudNH8sD63l+BZWZWThwYBXjzGY28sPMArXsOZV2KmVlmHBgFWDxrDICPMsxsQHNgFGBG41AmjhjMQ887MMxs4HJgFEAS508fxeObXqPDl9ea2QDlwCjQudNH8fqh47yw80DWpZiZZcKBUaBzp40E4LGXPLmSmQ1MDowCTRwxmPENg3hso2/gM7OByYFRIEmcO30Uj730mkevNbMByYHRC4umjeTVA0fZtNv3Y5jZwOPA6IWzJg8HYPWWPRlXYmbW/xwYvTBzTD1Dayp58uXXsy7FzKzfOTB6obJC/NbEBlZvcWCY2cDjwOilsyaPYN22fRw53p51KWZm/cqB0UsLJw3neHuwduuJs82amZW3NGfcu1nSTklruln/0byJldZIapc0Mlm3SdIzybqWtGo8FWdN6uz49mkpMxtY0jzCuAVY0t3KiPjHiFgYEQuBjwEPRkT+XXEXJ+ubU6yx18YMG8T4hkEODDMbcFILjIh4CCj0tuirgNvTqqWvLZw8nCdf9qW1ZjawZN6HIWkIuSOR7+YtDuB+SaskLcumsu6dOXE4rXsOs+fgsaxLMTPrN5kHBvBW4OcnnI66MCLOBi4DrpV0UXcbS1omqUVSy65d/TNfxfzxDQA8u80d32Y2cBRDYCzlhNNREbE1+boTuAtY1N3GEbE8IpojormxsTHVQjvNGz8MgDWv7O2XzzMzKwaZBoakBuDNwPfzlg2VVN/5HLgU6PJKq6yMGFrDhOGDfWmtmQ0oVWm9saTbgcXAaEmtwCeBaoCIuDFp9nbg/og4mLfpWOAuSZ31fSsi7k2rzlM1d/ww1mz1EYaZDRypBUZEXFVAm1vIXX6bv2wjsCCdqvrO/PEN/GjdDg4ebWNobWq70cysaBRDH0ZJmjd+GBGwzh3fZjZAODBO0fwJuSul3PFtZgOFA+MUjR1Wy6ihNe74NrMBw4FxiiQxb0IDaxwYZjZAODBOw7zxw3hhx36OtnmoczMrfw6M0zBv/DDaOoIXdhzIuhQzs9Q5ME7DnHG5O76f274/40rMzNLnwDgNU0cNpbaqgud8aa2ZDQAOjNNQWSHOGFvvIwwzGxAcGKdpdlM9z233EYaZlT8HxmmaPW4Yrx44xq79R7MuxcwsVQ6M0zSnqR6A9T4tZWZlzoFxmmYlgeHTUmZW7hwYp2lUXS1j6mtZt81HGGZW3hwYfWCWO77NbABwYPSBOeOG8cLOA7S1d2RdiplZalILDEk3S9opqcvpVSUtlrRX0urk8Ym8dUskrZe0QdL1adXYV2Y31XOsrYNNuw/23NjMrESleYRxC7CkhzYPR8TC5PFpAEmVwA3AZcBc4CpJc1Os87TNbsoNEeJ+DDMrZ6kFRkQ8BLx2CpsuAjZExMaIOAbcAVzRp8X1sRljhlJVIfdjmFlZy7oP43xJT0m6R9K8ZNkEYEtem9ZkWZckLZPUIqll165dadbardqqSmY01vGcjzDMrIxlGRhPAFMiYgHwr8D3kuXqom109yYRsTwimiOiubGxMYUyC5O7UsqBYWblK7PAiIh9EXEgeb4CqJY0mtwRxaS8phOBrRmU2Cuzx9XzyuuH2XfkeNalmJmlIrPAkNQkScnzRUktu4GVwExJ0yTVAEuBu7Oqs1Bzko5vn5Yys3JVldYbS7odWAyMltQKfBKoBoiIG4ErgQ9KagMOA0sjIoA2SdcB9wGVwM0RsTatOvvKrF+OKbWPRdNGZlyNmVnfSy0wIuKqHtZ/EfhiN+tWACvSqCst4xoGUT+oyv0YZla2sr5KqmxIYk7TMI9aa2Zly4HRh2Y11bN++35yZ9bMzMqLA6MPzWqqZ//RNl55/XDWpZiZ9TkHRh+aM86TKZlZ+XJg9KEzxnZOpuTAMLPy48DoQ/WDqpkwfLADw8zKkgOjj80ZV896D0JoZmXIgdHHZjXV8+Kugxxta8+6FDOzPuXA6GOzmobR3hG8uNOTKZlZeXFg9LE5nUOE7PBpKTMrLw6MPjZ19FBqKis8CKGZlR0HRh+rrqxgxpg6XyllZmXHgZGCOckQIWZm5cSBkYJZTfVs33eE1w8dy7oUM7M+48BIQefcGD4tZWblJLXAkHSzpJ2S1nSz/l2Snk4ej0hakLduk6RnJK2W1JJWjWmZMy43+55PS5lZOUnzCOMWYMlJ1r8EvDkizgT+Flh+wvqLI2JhRDSnVF9qxtTXMnxINc/5jm8zKyNpzrj3kKSpJ1n/SN7LR4GJadXS3yQxa2y9T0mZWVkplj6M9wP35L0O4H5JqyQtO9mGkpZJapHUsmvXrlSL7I3ZTfU8v30/HR2eTMnMykPmgSHpYnKB8Vd5iy+MiLOBy4BrJV3U3fYRsTwimiOiubGxMeVqCzd73DAOHmundY8nUzKz8pBpYEg6E7gJuCIidncuj4itydedwF3AomwqPHW/ulLK/RhmVh4yCwxJk4E7gfdExPN5y4dKqu98DlwKdHmlVTHrnEzJV0qZWblIrdNb0u3AYmC0pFbgk0A1QETcCHwCGAV8SRJAW3JF1FjgrmRZFfCtiLg3rTrTUldbxeSRQ9zxbWZlI82rpK7qYf0HgA90sXwjsOA3tyg9s5rqfUrKzMpG5p3e5Wx2Uz0vvXqQI8c9mZKZlT4HRopmNw2jI2DDzgNZl2JmdtocGCnymFJmVk4cGCmaOmoINVUVrHc/hpmVAQdGiqoqKzhjrCdTMrPy4MBI2ayxwxwYZlYWegwMSWMlfU3SPcnruZLen35p5WF2Uz279h/ltYOeTMnMSlshRxi3APcB45PXzwMfSaugcjN7nIcIMbPyUEhgjI6I7wAdABHRBvjGggL98kqpbT4tZWalrZDAOChpFLkhx5F0HrA31arKSGNdLSOH1nhMKTMreYUMDfLnwN3ADEk/BxqBK1Otqoz8cjKlHQ4MMyttPQZGRDwh6c3ALEDA+og4nnplZWT2uHrueHwLHR1BRYWyLsfM7JT0GBiSrj5h0dmSiIjbUqqp7Mxuqufw8XZefu0QU0cPzbocM7NTUsgpqTfmPR8EXAI8ATgwCjSraRiQGyLEgWFmpaqQU1Ifyn8tqQH4emoVlaEzxtYh5S6tXTK/KetyzMxOyanc6X0ImNlTI0k3S9opqcvZ8pTzBUkbJD0t6ey8dUskrU/WXX8KNRaVITVVTBk5xFdKmVlJK6QP4wckl9SSC5i5wHcKeO9bgC/S/amry8gFz0zgXODLwLmSKoEbgLcArcBKSXdHxLMFfGbRmtVU78Aws5JWSB/GZ/OetwGbI6K1p40i4iFJU0/S5ArgtogI4FFJwyWNA6YCG5KZ95B0R9K2pANjdtMw7n92B4ePtTO4pjLrcszMeq2QPowHU/rsCcCWvNetybKulp/b3ZtIWgYsA5g8eXLfV9lH5oyrJwKe37GfBZOGZ12OmVmvdduHIWm/pH1dPPZL6ouBkbq6ISFOsrxLEbE8IpojormxsbEPykrH3HENAKzd6jGlzKw0dXuEERH1KX92KzAp7/VEYCtQ083ykjZp5GDqB1WxZqtHVTGz0lTwVVKSxkia3Pnog8++G7g6uVrqPGBvRGwDVgIzJU2TVAMsTdqWNEnMH9/gIwwzK1mFXCX1B8A/kRvefCcwBVgHzOthu9uBxcBoSa3AJ4FqgIi4EVgBXA5sIHep7vuSdW2SriM3pHolcHNErD2F763ozBs/jNse3czx9g6qKz13lZmVlkKukvpb4DzgRxFxlqSLgat62igiTtomuTrq2m7WrSAXKGVl/oQGjrV18OKuA8xO7v42MysVhfyZezwidgMVkioi4qfAwpTrKkvzJ+RCYu0rPi1lZqWnkMB4XVId8DDwTUmfJ3c/hvXStNF1DK6udMe3mZWkQgLjIWA48GHgXuBF4K1pFlWuKivEnHH1PsIws5JUSGCIXAf0A0Ad8O3kFJWdgvkTGnh22z46Orq9tcTMrCj1GBgR8TcRMY9cB/V44EFJP0q9sjI1f3wDB462sfm1Q1mXYmbWK725tnMnsB3YDYxJp5zyN3d8ruN7zSvuxzCz0tJjYEj6oKQHgB8Do4H/GhFnpl1YuTpjbD3VlfINfGZWcgq5D2MK8JGIWJ12MQNBTVUFs5rqWesrpcysxBTSh3G9w6JvzR/fwJpX9pK7d9HMrDR4fIoMzJ/QwJ5Dx2ndczjrUszMCubAyMDCZD6Mp1pfz7gSM7PCOTAyMKupnpqqCp7a4sAws9LhwMhAdWUF88cP46kt7vg2s9LhwMjIgknDeeaVvbS1d2RdiplZQRwYGVk4aTiHj7fzws4DWZdiZlaQVAND0hJJ6yVtkHR9F+s/Kml18lgjqV3SyGTdJknPJOta0qwzCwsmJh3f7scwsxKRWmBIqgRuAC4D5gJXSZqb3yYi/jEiFkbEQuBjwIMR8Vpek4uT9c1p1ZmVKaOG0DC42ldKmVnJSPMIYxGwISI2RsQx4A7gipO0vwq4PcV6iookFkwazmp3fJtZiUgzMCYAW/JetybLfoOkIcAS4Lt5iwO4X9IqScu6+xBJyyS1SGrZtWtXH5TdfxZObGD99n0cOub5qMys+KUZGOpiWXdjYbwV+PkJp6MujIizyZ3SulbSRV1tGBHLI6I5IpobGxtPr+J+tmDScDoC1nhCJTMrAWkGRiswKe/1RGBrN22XcsLpqIjYmnzdCdxF7hRXWTnTHd9mVkLSDIyVwExJ0yTVkAuFu09sJKkBeDPw/bxlQyXVdz4HLgXWpFhrJhrra5k4YjBPbtmTdSlmZj0qZHjzUxIRbZKuIze9ayVwc0SslXRNsv7GpOnbgfsj4mDe5mOBuyR11vitiLg3rVqz1DxlBI+8uJuIIPl+zcyKUmqBARARK4AVJyy78YTXtwC3nLBsI7AgzdqKxTlTR/K91VvZ8tphJo8aknU5Zmbd8p3eGWueMgKAls2v9dDSzCxbDoyMnTG2nvpBVazc5H4MMytuDoyMVVaIsyePYJWPMMysyDkwisAbp47g+R0H2HvoeNalmJl1y4FRBM6ZMhKAVS/7KMPMipcDowgsnDScqgrR4n4MMytiDowiMLimknkTGhwYZlbUHBhFonnKCJ5qfZ2jbe1Zl2Jm1iUHRpE4d9pIjrZ1sPpljytlZsXJgVEkzp0+igrBIy/uzroUM7MuOTCKRMPgauZPaOAXGx0YZlacHBhF5Pzpo3jy5T0cPuZ+DDMrPg6MInL+jFEcbw+PK2VmRcmBUUTeOHUkVRVyP4aZFSUHRhEZWlvFgknD+YUDw8yKUKqBIWmJpPWSNki6vov1iyXtlbQ6eXyi0G3L1QUzRvF06+vsO+JxpcysuKQWGJIqgRuAy4C5wFWS5nbR9OGIWJg8Pt3LbcvO+TNG0RHw2Eb3Y5hZcUnzCGMRsCEiNkbEMeAO4Ip+2LaknTNlBENqKnnw+Z1Zl2Jm9mvSDIwJwJa8163JshOdL+kpSfdImtfLbctObVUlF8wYzQPrdxERWZdjZvZLaQaGulh24v+ATwBTImIB8K/A93qxba6htExSi6SWXbt2nXKxxWTxrEZa9xzmxV0Hsy7FzOyX0gyMVmBS3uuJwNb8BhGxLyIOJM9XANWSRheybd57LI+I5ohobmxs7Mv6M7N4Vu77eGC9T0uZWfFIMzBWAjMlTZNUAywF7s5vIKlJkpLni5J6dheybTmbOGIIM8fU8cD68jhiMrPyUJXWG0dEm6TrgPuASuDmiFgr6Zpk/Y3AlcAHJbUBh4GlkTtx3+W2adVajBbPauTWRzZz8GgbQ2tT+2cyMyuYyqljtbm5OVpaWrIuo088suFV3nnTY9x0dTO/O3ds1uWYWZmStCoimgtp6zu9i1Tz1JHU1Vbx4+d2ZF2KmRngwChaNVUVXDx7DPev3UF7R/kcBZpZ6XJgFLEl85rYffAYLZt817eZZc+BUcQWz2qktqqCe9duz7oUMzMHRjEbWlvFRWc0ct+a7b7r28wy58AockvmNbF17xGebt2bdSlmNsA5MIrcJXPGUFUhVqzZlnUpZjbAOTCK3PAhNfz2zNH8YPVWOny1lJllyIFRAt5+9kS27j3CYy/5aikzy44DowS8Zc5Y6mqruOvJ1qxLMbMBzIFRAgbXVLJkfhP3PLOdI8fbsy7HzAYoB0aJePtZE9h/tI0frfNQIWaWDQdGiThv+ijGNQzi2yu39NzYzCwFDowSUVkhlr5xMg+/8Cqbd3smPjPrfw6MEvKON06iskJ867GXsy7FzAYgB0YJaWoYxFvmjOU7LVs42ubObzPrX6kGhqQlktZL2iDp+i7Wv0vS08njEUkL8tZtkvSMpNWSymNWpD7w7vOmsOfQce55xgMSmln/Si0wJFUCNwCXAXOBqyTNPaHZS8CbI+JM4G+B5SesvzgiFhY6G9RAcMGMUUxvHMpXH97oAQnNrF+leYSxCNgQERsj4hhwB3BFfoOIeCQi9iQvHwUmplhPWaioEH960XTWbt3Hzza8mnU5ZjaApBkYE4D8a0Bbk2XdeT9wT97rAO6XtErSsu42krRMUoukll27dp1WwaXibWdNYOywWm588MWsSzGzASTNwFAXy7o8hyLpYnKB8Vd5iy+MiLPJndK6VtJFXW0bEcsjojkimhsbG0+35pJQW1XJn1w4jZ9v2M0zHvbczPpJmoHRCkzKez0R2HpiI0lnAjcBV0TE7s7lEbE1+boTuIvcKS5LvPPcyQwbVMW//Oj5rEsxswEizcBYCcyUNE1SDbAUuDu/gaTJwJ3AeyLi+bzlQyXVdz4HLgXWpFhryakfVM01i2fw4+d2stJzfptZP0gtMCKiDbgOuA9YB3wnItZKukbSNUmzTwCjgC+dcPnsWOBnkp4CHgf+f0Tcm1atpep9F0xjTH0t/3DPc75iysxSp3L6j6a5uTlaWgbWLRvffGwzH79rDcvfcw6XzmvKuhwzKzGSVhV664Lv9C5x/6V5EjPH1PHpHz7L4WO++9vM0uPAKHHVlRX83dvm07rnMP/6kxeyLsfMypgDowycO30U//nsiSx/aCPrt+/PuhwzK1MOjDLxPy6fTcPgaj58x5MemNDMUuHAKBOj6mr5zJVn8tz2/fzjveuzLsfMypADo4xcMmcs7zlvCjf97CXuW+vRbM2sbzkwyszHf28OCyY28GffXs26bfuyLsfMyogDo8wMqq5k+dXN1A+q4gO3trD19cNZl2RmZcKBUYbGDhvE1977RvYdOc47v/ooO/cdybokMysDDowyNX9CA7e8bxG79h9l6fJHeXn3oaxLMrMS58AoY+dMGcGtf7KI1w4d4+1f+jmrNu/peSMzs244MMpc89SR3PnBC6gbVMU7vvILvvzAi7R3lM/4YWbWfxwYA8D0xjq+f+2FvGXuWP7h3ue4avmjvoLKzHrNgTFADB9Sw5fedTafufJMnt+5n9/7wsN87M5n2Lz7YNalmVmJ8PDmA9Drh47xLz96gW899jJtHR0smd/EH50ziTfNHE11pf+GMBtIejO8eaqBIWkJ8HmgErgpIv7+hPVK1l8OHAL+OCKeKGTbrjgwemfHviP8359v4vbHX2bv4eOMHFrD4lmNXDBjNOdNH8mE4YPJ/ROZWbkqisCQVAk8D7yF3PzeK4GrIuLZvDaXAx8iFxjnAp+PiHML2bYrDoxTc6ytgwef38UPntrKzza8ymsHjwEwbFAVs8cN44yxdYwfPpimYYNoahjE6Lpa6mqrGFpbRV1tFZUVDhWzUtWbwKhKsY5FwIaI2JgUdQdwBZD/n/4VwG2RS61HJQ2XNA6YWsC21kdqqip4y9yxvGXuWDo6gue272fV5tdYt30/z23bx92rt7LvSFu32w+qrqC6soKqClFZ0flVVFeKCgl6yJOe4qanoxzHlQ10I4bU8J1rzk/9c9IMjAnAlrzXreSOInpqM6HAbQGQtAxYBjB58uTTq9ioqBBzxw9j7vhhv7b80LE2tu89wva9R9h98BgHj7ZxIHkcOtbO8fYO2juCto6gvT043pF73dMlvD0e3/bQIHp+B7OyN2xQdb98TpqB0dUffif+dnfXppBtcwsjlgPLIXdKqjcFWuGG1FQxvbGO6Y11WZdiZhlJMzBagUl5rycCWwtsU1PAtmZm1o/SvIZyJTBT0jRJNcBS4O4T2twNXK2c84C9EbGtwG3NzKwfpXaEERFtkq4D7iN3aezNEbFW0jXJ+huBFeSukNpA7rLa951s27RqNTOznvnGPTOzAaw3l9X6tl4zMyuIA8PMzAriwDAzs4I4MMzMrCBl1ektaRew+RQ3Hw282ofl9BXX1Tuuq3dcV++UY11TIqKxkIZlFRinQ1JLoVcK9CfX1Tuuq3dcV+8M9Lp8SsrMzAriwDAzs4I4MH5ledYFdMN19Y7r6h3X1TsDui73YZiZWUF8hGFmZgVxYJiZWUEGfGBIWiJpvaQNkq7P4PM3SXpG0mpJLcmykZL+XdILydcRee0/ltS6XtJ/6sM6bpa0U9KavGW9rkPSOcn3s0HSF9TT/KqnVtenJL2S7LPVydzw/V3XJEk/lbRO0lpJH06WZ7rPTlJXpvtM0iBJj0t6Kqnrb5LlWe+v7urK/Gcsec9KSU9K+mHyOtvfyYgYsA9yQ6e/CEwnN2nTU8Dcfq5hEzD6hGWfAa5Pnl8P/EPyfG5SYy0wLam9so/quAg4G1hzOnUAjwPnk5s18R7gshTq+hTwF1207c+6xgFnJ8/rgeeTz890n52krkz3WfIedcnzauAx4Lwi2F/d1ZX5z1jynn8OfAv4YTH8Tg7hTDyoAAAE5ElEQVT0I4xFwIaI2BgRx4A7gCsyrglyNdyaPL8VeFve8jsi4mhEvERuHpFFffGBEfEQ8Nrp1CFpHDAsIn4RuZ/U2/K26cu6utOfdW2LiCeS5/uBdeTmos90n52kru70V10REQeSl9XJI8h+f3VXV3f67WdM0kTg94CbTvj8zPbXQA+MCcCWvNetnPyXKw0B3C9plaRlybKxkZt5kOTrmGR5f9fb2zomJM/7o77rJD2dnLLqPCzPpC5JU4GzyP11WjT77IS6ION9lpxeWQ3sBP49Iopif3VTF2T/M/YvwF8CHXnLMt1fAz0wujqX19/XGV8YEWcDlwHXSrroJG2LoV7ovo7+qu/LwAxgIbAN+Kes6pJUB3wX+EhE7DtZ0/6srYu6Mt9nEdEeEQuBieT++p1/kuZZ15Xp/pL0+8DOiFhV6Cb9UddAD4xWYFLe64nA1v4sICK2Jl93AneRO8W0IzmUJPm6M2ne3/X2to7W5Hmq9UXEjuSXvAP4Kr86LdevdUmqJvef8jcj4s5kceb7rKu6imWfJbW8DjwALKEI9ldXdRXB/roQ+ANJm8idKv8dSd8g4/010ANjJTBT0jRJNcBS4O7++nBJQyXVdz4HLgXWJDW8N2n2XuD7yfO7gaWSaiVNA2aS69BKS6/qSA6R90s6L7kS4+q8bfpM5y9M4u3k9lm/1pW8z9eAdRHxubxVme6z7urKep9JapQ0PHk+GPhd4Dmy319d1pX1/oqIj0XExIiYSu7/pZ9ExLvJ+nfyVHvLy+UBXE7uSpIXgY/382dPJ3dlw1PA2s7PB0YBPwZeSL6OzNvm40mt6+mDqzDy3vd2cofex8n9VfL+U6kDaCb3y/Ui8EWS0QT6uK6vA88ATye/KOMyqOtN5A7tnwZWJ4/Ls95nJ6kr030GnAk8mXz+GuATp/qz3k91Zf4zlve+i/nVVVKZ7i8PDWJmZgUZ6KekzMysQA4MMzMriAPDzMwK4sAwM7OCODDMzKwgDgyzDEla3DkSqVmxc2CYmVlBHBhmBZD0buXmTVgt6SvJgHUHJP2TpCck/VhSY9J2oaRHk4Hr7uocuE7SGyT9SLm5F56QNCN5+zpJ/ybpOUnf7JyvQNLfS3o2eZ/PZvStm/2SA8OsB5LmAO8gN1DkQqAdeBcwFHgicoNHPgh8MtnkNuCvIuJMcncLdy7/JnBDRCwALiB3BzvkRpT9CLk5DaYDF0oaSW5IinnJ+/xdut+lWc8cGGY9uwQ4B1iZDIN9Cbn/2DuAbydtvgG8SVIDMDwiHkyW3wpclIwZNiEi7gKIiCMRcShp83hEtEZuoLvVwFRgH3AEuEnSHwKdbc0y48Aw65mAWyNiYfKYFRGf6qLdycbZOdm0mEfznrcDVRHRRm6E1O+Sm/Dm3l7WbNbnHBhmPfsxcKWkMfDLeZWnkPv9uTJp807gZxGxF9gj6beT5e8BHozcnBStkt6WvEetpCHdfWAyn0VDRKwgd7pqYRrfmFlvVGVdgFmxi4hnJf01uZkRK8iNnHstcBCYJ2kVsJdcPwfkhp2+MQmEjcD7kuXvAb4i6dPJe/zRST62Hvi+pEHkjk7+rI+/LbNe82i1ZqdI0oGIqMu6DrP+4lNSZmZWEB9hmJlZQXyEYWZmBXFgmJlZQRwYZmZWEAeGmZkVxIFhZmYF+Q9trXq8gf2FlAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(costs)\n",
    "plt.title('costs')\n",
    "plt.xlabel('epochs')\n",
    "plt.ylabel('value')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "951daa5e1959839fcb325fff331f52e72634f7a1be998f6081ed7f433b63f1b3"
  },
  "kernelspec": {
   "display_name": "Python 3.8.11 64-bit ('base': conda)",
   "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.8.11"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}