BMF_MCMC.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "BMF_MCMC.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyMB+Zi3dnpHgM3ynM+iPcG0",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "accelerator": "GPU"
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/ayush12gupta/Probabilistic-ML-ACA-/blob/master/BMF_MCMC.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "BL0YzfC8w8Hj",
        "colab_type": "code",
        "outputId": "c9dc8be7-ef87-4c02-d213-daee04931d17",
        "colab": {
          "resources": {
            "http://localhost:8080/nbextensions/google.colab/files.js": {
              "data": "Ly8gQ29weXJpZ2h0IDIwMTcgR29vZ2xlIExMQwovLwovLyBMaWNlbnNlZCB1bmRlciB0aGUgQXBhY2hlIExpY2Vuc2UsIFZlcnNpb24gMi4wICh0aGUgIkxpY2Vuc2UiKTsKLy8geW91IG1heSBub3QgdXNlIHRoaXMgZmlsZSBleGNlcHQgaW4gY29tcGxpYW5jZSB3aXRoIHRoZSBMaWNlbnNlLgovLyBZb3UgbWF5IG9idGFpbiBhIGNvcHkgb2YgdGhlIExpY2Vuc2UgYXQKLy8KLy8gICAgICBodHRwOi8vd3d3LmFwYWNoZS5vcmcvbGljZW5zZXMvTElDRU5TRS0yLjAKLy8KLy8gVW5sZXNzIHJlcXVpcmVkIGJ5IGFwcGxpY2FibGUgbGF3IG9yIGFncmVlZCB0byBpbiB3cml0aW5nLCBzb2Z0d2FyZQovLyBkaXN0cmlidXRlZCB1bmRlciB0aGUgTGljZW5zZSBpcyBkaXN0cmlidXRlZCBvbiBhbiAiQVMgSVMiIEJBU0lTLAovLyBXSVRIT1VUIFdBUlJBTlRJRVMgT1IgQ09ORElUSU9OUyBPRiBBTlkgS0lORCwgZWl0aGVyIGV4cHJlc3Mgb3IgaW1wbGllZC4KLy8gU2VlIHRoZSBMaWNlbnNlIGZvciB0aGUgc3BlY2lmaWMgbGFuZ3VhZ2UgZ292ZXJuaW5nIHBlcm1pc3Npb25zIGFuZAovLyBsaW1pdGF0aW9ucyB1bmRlciB0aGUgTGljZW5zZS4KCi8qKgogKiBAZmlsZW92ZXJ2aWV3IEhlbHBlcnMgZm9yIGdvb2dsZS5jb2xhYiBQeXRob24gbW9kdWxlLgogKi8KKGZ1bmN0aW9uKHNjb3BlKSB7CmZ1bmN0aW9uIHNwYW4odGV4dCwgc3R5bGVBdHRyaWJ1dGVzID0ge30pIHsKICBjb25zdCBlbGVtZW50ID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnc3BhbicpOwogIGVsZW1lbnQudGV4dENvbnRlbnQgPSB0ZXh0OwogIGZvciAoY29uc3Qga2V5IG9mIE9iamVjdC5rZXlzKHN0eWxlQXR0cmlidXRlcykpIHsKICAgIGVsZW1lbnQuc3R5bGVba2V5XSA9IHN0eWxlQXR0cmlidXRlc1trZXldOwogIH0KICByZXR1cm4gZWxlbWVudDsKfQoKLy8gTWF4IG51bWJlciBvZiBieXRlcyB3aGljaCB3aWxsIGJlIHVwbG9hZGVkIGF0IGEgdGltZS4KY29uc3QgTUFYX1BBWUxPQURfU0laRSA9IDEwMCAqIDEwMjQ7Ci8vIE1heCBhbW91bnQgb2YgdGltZSB0byBibG9jayB3YWl0aW5nIGZvciB0aGUgdXNlci4KY29uc3QgRklMRV9DSEFOR0VfVElNRU9VVF9NUyA9IDMwICogMTAwMDsKCmZ1bmN0aW9uIF91cGxvYWRGaWxlcyhpbnB1dElkLCBvdXRwdXRJZCkgewogIGNvbnN0IHN0ZXBzID0gdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKTsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIC8vIENhY2hlIHN0ZXBzIG9uIHRoZSBvdXRwdXRFbGVtZW50IHRvIG1ha2UgaXQgYXZhaWxhYmxlIGZvciB0aGUgbmV4dCBjYWxsCiAgLy8gdG8gdXBsb2FkRmlsZXNDb250aW51ZSBmcm9tIFB5dGhvbi4KICBvdXRwdXRFbGVtZW50LnN0ZXBzID0gc3RlcHM7CgogIHJldHVybiBfdXBsb2FkRmlsZXNDb250aW51ZShvdXRwdXRJZCk7Cn0KCi8vIFRoaXMgaXMgcm91Z2hseSBhbiBhc3luYyBnZW5lcmF0b3IgKG5vdCBzdXBwb3J0ZWQgaW4gdGhlIGJyb3dzZXIgeWV0KSwKLy8gd2hlcmUgdGhlcmUgYXJlIG11bHRpcGxlIGFzeW5jaHJvbm91cyBzdGVwcyBhbmQgdGhlIFB5dGhvbiBzaWRlIGlzIGdvaW5nCi8vIHRvIHBvbGwgZm9yIGNvbXBsZXRpb24gb2YgZWFjaCBzdGVwLgovLyBUaGlzIHVzZXMgYSBQcm9taXNlIHRvIGJsb2NrIHRoZSBweXRob24gc2lkZSBvbiBjb21wbGV0aW9uIG9mIGVhY2ggc3RlcCwKLy8gdGhlbiBwYXNzZXMgdGhlIHJlc3VsdCBvZiB0aGUgcHJldmlvdXMgc3RlcCBhcyB0aGUgaW5wdXQgdG8gdGhlIG5leHQgc3RlcC4KZnVuY3Rpb24gX3VwbG9hZEZpbGVzQ29udGludWUob3V0cHV0SWQpIHsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIGNvbnN0IHN0ZXBzID0gb3V0cHV0RWxlbWVudC5zdGVwczsKCiAgY29uc3QgbmV4dCA9IHN0ZXBzLm5leHQob3V0cHV0RWxlbWVudC5sYXN0UHJvbWlzZVZhbHVlKTsKICByZXR1cm4gUHJvbWlzZS5yZXNvbHZlKG5leHQudmFsdWUucHJvbWlzZSkudGhlbigodmFsdWUpID0+IHsKICAgIC8vIENhY2hlIHRoZSBsYXN0IHByb21pc2UgdmFsdWUgdG8gbWFrZSBpdCBhdmFpbGFibGUgdG8gdGhlIG5leHQKICAgIC8vIHN0ZXAgb2YgdGhlIGdlbmVyYXRvci4KICAgIG91dHB1dEVsZW1lbnQubGFzdFByb21pc2VWYWx1ZSA9IHZhbHVlOwogICAgcmV0dXJuIG5leHQudmFsdWUucmVzcG9uc2U7CiAgfSk7Cn0KCi8qKgogKiBHZW5lcmF0b3IgZnVuY3Rpb24gd2hpY2ggaXMgY2FsbGVkIGJldHdlZW4gZWFjaCBhc3luYyBzdGVwIG9mIHRoZSB1cGxvYWQKICogcHJvY2Vzcy4KICogQHBhcmFtIHtzdHJpbmd9IGlucHV0SWQgRWxlbWVudCBJRCBvZiB0aGUgaW5wdXQgZmlsZSBwaWNrZXIgZWxlbWVudC4KICogQHBhcmFtIHtzdHJpbmd9IG91dHB1dElkIEVsZW1lbnQgSUQgb2YgdGhlIG91dHB1dCBkaXNwbGF5LgogKiBAcmV0dXJuIHshSXRlcmFibGU8IU9iamVjdD59IEl0ZXJhYmxlIG9mIG5leHQgc3RlcHMuCiAqLwpmdW5jdGlvbiogdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKSB7CiAgY29uc3QgaW5wdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQoaW5wdXRJZCk7CiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gZmFsc2U7CgogIGNvbnN0IG91dHB1dEVsZW1lbnQgPSBkb2N1bWVudC5nZXRFbGVtZW50QnlJZChvdXRwdXRJZCk7CiAgb3V0cHV0RWxlbWVudC5pbm5lckhUTUwgPSAnJzsKCiAgY29uc3QgcGlja2VkUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBpbnB1dEVsZW1lbnQuYWRkRXZlbnRMaXN0ZW5lcignY2hhbmdlJywgKGUpID0+IHsKICAgICAgcmVzb2x2ZShlLnRhcmdldC5maWxlcyk7CiAgICB9KTsKICB9KTsKCiAgY29uc3QgY2FuY2VsID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnYnV0dG9uJyk7CiAgaW5wdXRFbGVtZW50LnBhcmVudEVsZW1lbnQuYXBwZW5kQ2hpbGQoY2FuY2VsKTsKICBjYW5jZWwudGV4dENvbnRlbnQgPSAnQ2FuY2VsIHVwbG9hZCc7CiAgY29uc3QgY2FuY2VsUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBjYW5jZWwub25jbGljayA9ICgpID0+IHsKICAgICAgcmVzb2x2ZShudWxsKTsKICAgIH07CiAgfSk7CgogIC8vIENhbmNlbCB1cGxvYWQgaWYgdXNlciBoYXNuJ3QgcGlja2VkIGFueXRoaW5nIGluIHRpbWVvdXQuCiAgY29uc3QgdGltZW91dFByb21pc2UgPSBuZXcgUHJvbWlzZSgocmVzb2x2ZSkgPT4gewogICAgc2V0VGltZW91dCgoKSA9PiB7CiAgICAgIHJlc29sdmUobnVsbCk7CiAgICB9LCBGSUxFX0NIQU5HRV9USU1FT1VUX01TKTsKICB9KTsKCiAgLy8gV2FpdCBmb3IgdGhlIHVzZXIgdG8gcGljayB0aGUgZmlsZXMuCiAgY29uc3QgZmlsZXMgPSB5aWVsZCB7CiAgICBwcm9taXNlOiBQcm9taXNlLnJhY2UoW3BpY2tlZFByb21pc2UsIHRpbWVvdXRQcm9taXNlLCBjYW5jZWxQcm9taXNlXSksCiAgICByZXNwb25zZTogewogICAgICBhY3Rpb246ICdzdGFydGluZycsCiAgICB9CiAgfTsKCiAgaWYgKCFmaWxlcykgewogICAgcmV0dXJuIHsKICAgICAgcmVzcG9uc2U6IHsKICAgICAgICBhY3Rpb246ICdjb21wbGV0ZScsCiAgICAgIH0KICAgIH07CiAgfQoKICBjYW5jZWwucmVtb3ZlKCk7CgogIC8vIERpc2FibGUgdGhlIGlucHV0IGVsZW1lbnQgc2luY2UgZnVydGhlciBwaWNrcyBhcmUgbm90IGFsbG93ZWQuCiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gdHJ1ZTsKCiAgZm9yIChjb25zdCBmaWxlIG9mIGZpbGVzKSB7CiAgICBjb25zdCBsaSA9IGRvY3VtZW50LmNyZWF0ZUVsZW1lbnQoJ2xpJyk7CiAgICBsaS5hcHBlbmQoc3BhbihmaWxlLm5hbWUsIHtmb250V2VpZ2h0OiAnYm9sZCd9KSk7CiAgICBsaS5hcHBlbmQoc3BhbigKICAgICAgICBgKCR7ZmlsZS50eXBlIHx8ICduL2EnfSkgLSAke2ZpbGUuc2l6ZX0gYnl0ZXMsIGAgKwogICAgICAgIGBsYXN0IG1vZGlmaWVkOiAkewogICAgICAgICAgICBmaWxlLmxhc3RNb2RpZmllZERhdGUgPyBmaWxlLmxhc3RNb2RpZmllZERhdGUudG9Mb2NhbGVEYXRlU3RyaW5nKCkgOgogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAnbi9hJ30gLSBgKSk7CiAgICBjb25zdCBwZXJjZW50ID0gc3BhbignMCUgZG9uZScpOwogICAgbGkuYXBwZW5kQ2hpbGQocGVyY2VudCk7CgogICAgb3V0cHV0RWxlbWVudC5hcHBlbmRDaGlsZChsaSk7CgogICAgY29uc3QgZmlsZURhdGFQcm9taXNlID0gbmV3IFByb21pc2UoKHJlc29sdmUpID0+IHsKICAgICAgY29uc3QgcmVhZGVyID0gbmV3IEZpbGVSZWFkZXIoKTsKICAgICAgcmVhZGVyLm9ubG9hZCA9IChlKSA9PiB7CiAgICAgICAgcmVzb2x2ZShlLnRhcmdldC5yZXN1bHQpOwogICAgICB9OwogICAgICByZWFkZXIucmVhZEFzQXJyYXlCdWZmZXIoZmlsZSk7CiAgICB9KTsKICAgIC8vIFdhaXQgZm9yIHRoZSBkYXRhIHRvIGJlIHJlYWR5LgogICAgbGV0IGZpbGVEYXRhID0geWllbGQgewogICAgICBwcm9taXNlOiBmaWxlRGF0YVByb21pc2UsCiAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgYWN0aW9uOiAnY29udGludWUnLAogICAgICB9CiAgICB9OwoKICAgIC8vIFVzZSBhIGNodW5rZWQgc2VuZGluZyB0byBhdm9pZCBtZXNzYWdlIHNpemUgbGltaXRzLiBTZWUgYi82MjExNTY2MC4KICAgIGxldCBwb3NpdGlvbiA9IDA7CiAgICB3aGlsZSAocG9zaXRpb24gPCBmaWxlRGF0YS5ieXRlTGVuZ3RoKSB7CiAgICAgIGNvbnN0IGxlbmd0aCA9IE1hdGgubWluKGZpbGVEYXRhLmJ5dGVMZW5ndGggLSBwb3NpdGlvbiwgTUFYX1BBWUxPQURfU0laRSk7CiAgICAgIGNvbnN0IGNodW5rID0gbmV3IFVpbnQ4QXJyYXkoZmlsZURhdGEsIHBvc2l0aW9uLCBsZW5ndGgpOwogICAgICBwb3NpdGlvbiArPSBsZW5ndGg7CgogICAgICBjb25zdCBiYXNlNjQgPSBidG9hKFN0cmluZy5mcm9tQ2hhckNvZGUuYXBwbHkobnVsbCwgY2h1bmspKTsKICAgICAgeWllbGQgewogICAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgICBhY3Rpb246ICdhcHBlbmQnLAogICAgICAgICAgZmlsZTogZmlsZS5uYW1lLAogICAgICAgICAgZGF0YTogYmFzZTY0LAogICAgICAgIH0sCiAgICAgIH07CiAgICAgIHBlcmNlbnQudGV4dENvbnRlbnQgPQogICAgICAgICAgYCR7TWF0aC5yb3VuZCgocG9zaXRpb24gLyBmaWxlRGF0YS5ieXRlTGVuZ3RoKSAqIDEwMCl9JSBkb25lYDsKICAgIH0KICB9CgogIC8vIEFsbCBkb25lLgogIHlpZWxkIHsKICAgIHJlc3BvbnNlOiB7CiAgICAgIGFjdGlvbjogJ2NvbXBsZXRlJywKICAgIH0KICB9Owp9CgpzY29wZS5nb29nbGUgPSBzY29wZS5nb29nbGUgfHwge307CnNjb3BlLmdvb2dsZS5jb2xhYiA9IHNjb3BlLmdvb2dsZS5jb2xhYiB8fCB7fTsKc2NvcGUuZ29vZ2xlLmNvbGFiLl9maWxlcyA9IHsKICBfdXBsb2FkRmlsZXMsCiAgX3VwbG9hZEZpbGVzQ29udGludWUsCn07Cn0pKHNlbGYpOwo=",
              "ok": true,
              "headers": [
                [
                  "content-type",
                  "application/javascript"
                ]
              ],
              "status": 200,
              "status_text": ""
            }
          },
          "base_uri": "https://localhost:8080/",
          "height": 74
        }
      },
      "source": [
        "# Importing Train File\n",
        "from google.colab import files\n",
        "uploaded = files.upload()\n"
      ],
      "execution_count": 6,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/html": [
              "\n",
              "     <input type=\"file\" id=\"files-f06169c8-0676-4e35-8038-5db2ba02799d\" name=\"files[]\" multiple disabled />\n",
              "     <output id=\"result-f06169c8-0676-4e35-8038-5db2ba02799d\">\n",
              "      Upload widget is only available when the cell has been executed in the\n",
              "      current browser session. Please rerun this cell to enable.\n",
              "      </output>\n",
              "      <script src=\"/nbextensions/google.colab/files.js\"></script> "
            ],
            "text/plain": [
              "<IPython.core.display.HTML object>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "stream",
          "text": [
            "Saving u1.base to u1.base\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "fxtzg5g_xJH4",
        "colab_type": "code",
        "outputId": "f0805819-41eb-4dc7-dd3e-28b4b372b7c1",
        "colab": {
          "resources": {
            "http://localhost:8080/nbextensions/google.colab/files.js": {
              "data": "Ly8gQ29weXJpZ2h0IDIwMTcgR29vZ2xlIExMQwovLwovLyBMaWNlbnNlZCB1bmRlciB0aGUgQXBhY2hlIExpY2Vuc2UsIFZlcnNpb24gMi4wICh0aGUgIkxpY2Vuc2UiKTsKLy8geW91IG1heSBub3QgdXNlIHRoaXMgZmlsZSBleGNlcHQgaW4gY29tcGxpYW5jZSB3aXRoIHRoZSBMaWNlbnNlLgovLyBZb3UgbWF5IG9idGFpbiBhIGNvcHkgb2YgdGhlIExpY2Vuc2UgYXQKLy8KLy8gICAgICBodHRwOi8vd3d3LmFwYWNoZS5vcmcvbGljZW5zZXMvTElDRU5TRS0yLjAKLy8KLy8gVW5sZXNzIHJlcXVpcmVkIGJ5IGFwcGxpY2FibGUgbGF3IG9yIGFncmVlZCB0byBpbiB3cml0aW5nLCBzb2Z0d2FyZQovLyBkaXN0cmlidXRlZCB1bmRlciB0aGUgTGljZW5zZSBpcyBkaXN0cmlidXRlZCBvbiBhbiAiQVMgSVMiIEJBU0lTLAovLyBXSVRIT1VUIFdBUlJBTlRJRVMgT1IgQ09ORElUSU9OUyBPRiBBTlkgS0lORCwgZWl0aGVyIGV4cHJlc3Mgb3IgaW1wbGllZC4KLy8gU2VlIHRoZSBMaWNlbnNlIGZvciB0aGUgc3BlY2lmaWMgbGFuZ3VhZ2UgZ292ZXJuaW5nIHBlcm1pc3Npb25zIGFuZAovLyBsaW1pdGF0aW9ucyB1bmRlciB0aGUgTGljZW5zZS4KCi8qKgogKiBAZmlsZW92ZXJ2aWV3IEhlbHBlcnMgZm9yIGdvb2dsZS5jb2xhYiBQeXRob24gbW9kdWxlLgogKi8KKGZ1bmN0aW9uKHNjb3BlKSB7CmZ1bmN0aW9uIHNwYW4odGV4dCwgc3R5bGVBdHRyaWJ1dGVzID0ge30pIHsKICBjb25zdCBlbGVtZW50ID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnc3BhbicpOwogIGVsZW1lbnQudGV4dENvbnRlbnQgPSB0ZXh0OwogIGZvciAoY29uc3Qga2V5IG9mIE9iamVjdC5rZXlzKHN0eWxlQXR0cmlidXRlcykpIHsKICAgIGVsZW1lbnQuc3R5bGVba2V5XSA9IHN0eWxlQXR0cmlidXRlc1trZXldOwogIH0KICByZXR1cm4gZWxlbWVudDsKfQoKLy8gTWF4IG51bWJlciBvZiBieXRlcyB3aGljaCB3aWxsIGJlIHVwbG9hZGVkIGF0IGEgdGltZS4KY29uc3QgTUFYX1BBWUxPQURfU0laRSA9IDEwMCAqIDEwMjQ7Ci8vIE1heCBhbW91bnQgb2YgdGltZSB0byBibG9jayB3YWl0aW5nIGZvciB0aGUgdXNlci4KY29uc3QgRklMRV9DSEFOR0VfVElNRU9VVF9NUyA9IDMwICogMTAwMDsKCmZ1bmN0aW9uIF91cGxvYWRGaWxlcyhpbnB1dElkLCBvdXRwdXRJZCkgewogIGNvbnN0IHN0ZXBzID0gdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKTsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIC8vIENhY2hlIHN0ZXBzIG9uIHRoZSBvdXRwdXRFbGVtZW50IHRvIG1ha2UgaXQgYXZhaWxhYmxlIGZvciB0aGUgbmV4dCBjYWxsCiAgLy8gdG8gdXBsb2FkRmlsZXNDb250aW51ZSBmcm9tIFB5dGhvbi4KICBvdXRwdXRFbGVtZW50LnN0ZXBzID0gc3RlcHM7CgogIHJldHVybiBfdXBsb2FkRmlsZXNDb250aW51ZShvdXRwdXRJZCk7Cn0KCi8vIFRoaXMgaXMgcm91Z2hseSBhbiBhc3luYyBnZW5lcmF0b3IgKG5vdCBzdXBwb3J0ZWQgaW4gdGhlIGJyb3dzZXIgeWV0KSwKLy8gd2hlcmUgdGhlcmUgYXJlIG11bHRpcGxlIGFzeW5jaHJvbm91cyBzdGVwcyBhbmQgdGhlIFB5dGhvbiBzaWRlIGlzIGdvaW5nCi8vIHRvIHBvbGwgZm9yIGNvbXBsZXRpb24gb2YgZWFjaCBzdGVwLgovLyBUaGlzIHVzZXMgYSBQcm9taXNlIHRvIGJsb2NrIHRoZSBweXRob24gc2lkZSBvbiBjb21wbGV0aW9uIG9mIGVhY2ggc3RlcCwKLy8gdGhlbiBwYXNzZXMgdGhlIHJlc3VsdCBvZiB0aGUgcHJldmlvdXMgc3RlcCBhcyB0aGUgaW5wdXQgdG8gdGhlIG5leHQgc3RlcC4KZnVuY3Rpb24gX3VwbG9hZEZpbGVzQ29udGludWUob3V0cHV0SWQpIHsKICBjb25zdCBvdXRwdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQob3V0cHV0SWQpOwogIGNvbnN0IHN0ZXBzID0gb3V0cHV0RWxlbWVudC5zdGVwczsKCiAgY29uc3QgbmV4dCA9IHN0ZXBzLm5leHQob3V0cHV0RWxlbWVudC5sYXN0UHJvbWlzZVZhbHVlKTsKICByZXR1cm4gUHJvbWlzZS5yZXNvbHZlKG5leHQudmFsdWUucHJvbWlzZSkudGhlbigodmFsdWUpID0+IHsKICAgIC8vIENhY2hlIHRoZSBsYXN0IHByb21pc2UgdmFsdWUgdG8gbWFrZSBpdCBhdmFpbGFibGUgdG8gdGhlIG5leHQKICAgIC8vIHN0ZXAgb2YgdGhlIGdlbmVyYXRvci4KICAgIG91dHB1dEVsZW1lbnQubGFzdFByb21pc2VWYWx1ZSA9IHZhbHVlOwogICAgcmV0dXJuIG5leHQudmFsdWUucmVzcG9uc2U7CiAgfSk7Cn0KCi8qKgogKiBHZW5lcmF0b3IgZnVuY3Rpb24gd2hpY2ggaXMgY2FsbGVkIGJldHdlZW4gZWFjaCBhc3luYyBzdGVwIG9mIHRoZSB1cGxvYWQKICogcHJvY2Vzcy4KICogQHBhcmFtIHtzdHJpbmd9IGlucHV0SWQgRWxlbWVudCBJRCBvZiB0aGUgaW5wdXQgZmlsZSBwaWNrZXIgZWxlbWVudC4KICogQHBhcmFtIHtzdHJpbmd9IG91dHB1dElkIEVsZW1lbnQgSUQgb2YgdGhlIG91dHB1dCBkaXNwbGF5LgogKiBAcmV0dXJuIHshSXRlcmFibGU8IU9iamVjdD59IEl0ZXJhYmxlIG9mIG5leHQgc3RlcHMuCiAqLwpmdW5jdGlvbiogdXBsb2FkRmlsZXNTdGVwKGlucHV0SWQsIG91dHB1dElkKSB7CiAgY29uc3QgaW5wdXRFbGVtZW50ID0gZG9jdW1lbnQuZ2V0RWxlbWVudEJ5SWQoaW5wdXRJZCk7CiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gZmFsc2U7CgogIGNvbnN0IG91dHB1dEVsZW1lbnQgPSBkb2N1bWVudC5nZXRFbGVtZW50QnlJZChvdXRwdXRJZCk7CiAgb3V0cHV0RWxlbWVudC5pbm5lckhUTUwgPSAnJzsKCiAgY29uc3QgcGlja2VkUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBpbnB1dEVsZW1lbnQuYWRkRXZlbnRMaXN0ZW5lcignY2hhbmdlJywgKGUpID0+IHsKICAgICAgcmVzb2x2ZShlLnRhcmdldC5maWxlcyk7CiAgICB9KTsKICB9KTsKCiAgY29uc3QgY2FuY2VsID0gZG9jdW1lbnQuY3JlYXRlRWxlbWVudCgnYnV0dG9uJyk7CiAgaW5wdXRFbGVtZW50LnBhcmVudEVsZW1lbnQuYXBwZW5kQ2hpbGQoY2FuY2VsKTsKICBjYW5jZWwudGV4dENvbnRlbnQgPSAnQ2FuY2VsIHVwbG9hZCc7CiAgY29uc3QgY2FuY2VsUHJvbWlzZSA9IG5ldyBQcm9taXNlKChyZXNvbHZlKSA9PiB7CiAgICBjYW5jZWwub25jbGljayA9ICgpID0+IHsKICAgICAgcmVzb2x2ZShudWxsKTsKICAgIH07CiAgfSk7CgogIC8vIENhbmNlbCB1cGxvYWQgaWYgdXNlciBoYXNuJ3QgcGlja2VkIGFueXRoaW5nIGluIHRpbWVvdXQuCiAgY29uc3QgdGltZW91dFByb21pc2UgPSBuZXcgUHJvbWlzZSgocmVzb2x2ZSkgPT4gewogICAgc2V0VGltZW91dCgoKSA9PiB7CiAgICAgIHJlc29sdmUobnVsbCk7CiAgICB9LCBGSUxFX0NIQU5HRV9USU1FT1VUX01TKTsKICB9KTsKCiAgLy8gV2FpdCBmb3IgdGhlIHVzZXIgdG8gcGljayB0aGUgZmlsZXMuCiAgY29uc3QgZmlsZXMgPSB5aWVsZCB7CiAgICBwcm9taXNlOiBQcm9taXNlLnJhY2UoW3BpY2tlZFByb21pc2UsIHRpbWVvdXRQcm9taXNlLCBjYW5jZWxQcm9taXNlXSksCiAgICByZXNwb25zZTogewogICAgICBhY3Rpb246ICdzdGFydGluZycsCiAgICB9CiAgfTsKCiAgaWYgKCFmaWxlcykgewogICAgcmV0dXJuIHsKICAgICAgcmVzcG9uc2U6IHsKICAgICAgICBhY3Rpb246ICdjb21wbGV0ZScsCiAgICAgIH0KICAgIH07CiAgfQoKICBjYW5jZWwucmVtb3ZlKCk7CgogIC8vIERpc2FibGUgdGhlIGlucHV0IGVsZW1lbnQgc2luY2UgZnVydGhlciBwaWNrcyBhcmUgbm90IGFsbG93ZWQuCiAgaW5wdXRFbGVtZW50LmRpc2FibGVkID0gdHJ1ZTsKCiAgZm9yIChjb25zdCBmaWxlIG9mIGZpbGVzKSB7CiAgICBjb25zdCBsaSA9IGRvY3VtZW50LmNyZWF0ZUVsZW1lbnQoJ2xpJyk7CiAgICBsaS5hcHBlbmQoc3BhbihmaWxlLm5hbWUsIHtmb250V2VpZ2h0OiAnYm9sZCd9KSk7CiAgICBsaS5hcHBlbmQoc3BhbigKICAgICAgICBgKCR7ZmlsZS50eXBlIHx8ICduL2EnfSkgLSAke2ZpbGUuc2l6ZX0gYnl0ZXMsIGAgKwogICAgICAgIGBsYXN0IG1vZGlmaWVkOiAkewogICAgICAgICAgICBmaWxlLmxhc3RNb2RpZmllZERhdGUgPyBmaWxlLmxhc3RNb2RpZmllZERhdGUudG9Mb2NhbGVEYXRlU3RyaW5nKCkgOgogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAnbi9hJ30gLSBgKSk7CiAgICBjb25zdCBwZXJjZW50ID0gc3BhbignMCUgZG9uZScpOwogICAgbGkuYXBwZW5kQ2hpbGQocGVyY2VudCk7CgogICAgb3V0cHV0RWxlbWVudC5hcHBlbmRDaGlsZChsaSk7CgogICAgY29uc3QgZmlsZURhdGFQcm9taXNlID0gbmV3IFByb21pc2UoKHJlc29sdmUpID0+IHsKICAgICAgY29uc3QgcmVhZGVyID0gbmV3IEZpbGVSZWFkZXIoKTsKICAgICAgcmVhZGVyLm9ubG9hZCA9IChlKSA9PiB7CiAgICAgICAgcmVzb2x2ZShlLnRhcmdldC5yZXN1bHQpOwogICAgICB9OwogICAgICByZWFkZXIucmVhZEFzQXJyYXlCdWZmZXIoZmlsZSk7CiAgICB9KTsKICAgIC8vIFdhaXQgZm9yIHRoZSBkYXRhIHRvIGJlIHJlYWR5LgogICAgbGV0IGZpbGVEYXRhID0geWllbGQgewogICAgICBwcm9taXNlOiBmaWxlRGF0YVByb21pc2UsCiAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgYWN0aW9uOiAnY29udGludWUnLAogICAgICB9CiAgICB9OwoKICAgIC8vIFVzZSBhIGNodW5rZWQgc2VuZGluZyB0byBhdm9pZCBtZXNzYWdlIHNpemUgbGltaXRzLiBTZWUgYi82MjExNTY2MC4KICAgIGxldCBwb3NpdGlvbiA9IDA7CiAgICB3aGlsZSAocG9zaXRpb24gPCBmaWxlRGF0YS5ieXRlTGVuZ3RoKSB7CiAgICAgIGNvbnN0IGxlbmd0aCA9IE1hdGgubWluKGZpbGVEYXRhLmJ5dGVMZW5ndGggLSBwb3NpdGlvbiwgTUFYX1BBWUxPQURfU0laRSk7CiAgICAgIGNvbnN0IGNodW5rID0gbmV3IFVpbnQ4QXJyYXkoZmlsZURhdGEsIHBvc2l0aW9uLCBsZW5ndGgpOwogICAgICBwb3NpdGlvbiArPSBsZW5ndGg7CgogICAgICBjb25zdCBiYXNlNjQgPSBidG9hKFN0cmluZy5mcm9tQ2hhckNvZGUuYXBwbHkobnVsbCwgY2h1bmspKTsKICAgICAgeWllbGQgewogICAgICAgIHJlc3BvbnNlOiB7CiAgICAgICAgICBhY3Rpb246ICdhcHBlbmQnLAogICAgICAgICAgZmlsZTogZmlsZS5uYW1lLAogICAgICAgICAgZGF0YTogYmFzZTY0LAogICAgICAgIH0sCiAgICAgIH07CiAgICAgIHBlcmNlbnQudGV4dENvbnRlbnQgPQogICAgICAgICAgYCR7TWF0aC5yb3VuZCgocG9zaXRpb24gLyBmaWxlRGF0YS5ieXRlTGVuZ3RoKSAqIDEwMCl9JSBkb25lYDsKICAgIH0KICB9CgogIC8vIEFsbCBkb25lLgogIHlpZWxkIHsKICAgIHJlc3BvbnNlOiB7CiAgICAgIGFjdGlvbjogJ2NvbXBsZXRlJywKICAgIH0KICB9Owp9CgpzY29wZS5nb29nbGUgPSBzY29wZS5nb29nbGUgfHwge307CnNjb3BlLmdvb2dsZS5jb2xhYiA9IHNjb3BlLmdvb2dsZS5jb2xhYiB8fCB7fTsKc2NvcGUuZ29vZ2xlLmNvbGFiLl9maWxlcyA9IHsKICBfdXBsb2FkRmlsZXMsCiAgX3VwbG9hZEZpbGVzQ29udGludWUsCn07Cn0pKHNlbGYpOwo=",
              "ok": true,
              "headers": [
                [
                  "content-type",
                  "application/javascript"
                ]
              ],
              "status": 200,
              "status_text": ""
            }
          },
          "base_uri": "https://localhost:8080/",
          "height": 74
        }
      },
      "source": [
        "# Importing Test File\n",
        "from google.colab import files\n",
        "uploaded = files.upload()"
      ],
      "execution_count": 8,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/html": [
              "\n",
              "     <input type=\"file\" id=\"files-7306043d-7011-4d30-a35a-5e51afd74909\" name=\"files[]\" multiple disabled />\n",
              "     <output id=\"result-7306043d-7011-4d30-a35a-5e51afd74909\">\n",
              "      Upload widget is only available when the cell has been executed in the\n",
              "      current browser session. Please rerun this cell to enable.\n",
              "      </output>\n",
              "      <script src=\"/nbextensions/google.colab/files.js\"></script> "
            ],
            "text/plain": [
              "<IPython.core.display.HTML object>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "stream",
          "text": [
            "Saving u1.test to u1.test\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ahr2kt4o0u6l",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "from scipy.stats import wishart"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "jPBTCWt42w4B",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "N = 943 # No. of users\n",
        "M = 1682 # No. of movies  \n",
        "D = 50 "
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "RiDaiyHn1k_5",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def matrix_generator(N,M,file):\n",
        "\n",
        "    I = np.zeros((N,M))\n",
        "    R = np.zeros((N,M))\n",
        "    spt = \"\t\"\n",
        "    op = open(file,\"r\")\n",
        "    for line in op:\n",
        "        (usr,movie,rating,_) = line.split(spt)\n",
        "        #check !!!!!!!!\n",
        "        R[int(usr) - 1, int(movie) - 1] = int(rating)\n",
        "        I[int(usr) - 1, int(movie) - 1] = 1.0\n",
        "    return R, I"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "mPeQGan61EQ5",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "R, I = matrix_generator(N,M,\"u1.base\")\n",
        "RTest, ITest = matrix_generator(N,M,\"u1.test\")"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "RtTZKTZx_6ur",
        "colab_type": "text"
      },
      "source": [
        "**Hyperparameter Declaration**"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "8WEPOVoq_zYG",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "w0 = np.eye(D)\n",
        "mu0 = np.zeros(D)\n",
        "v0 = D\n",
        "alpha = 2.0\n",
        "beta0 = 1.0\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Pn8vmINiCmDW",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "lambv = wishart.rvs(size = 1,df = vv0, scale = wv0)\n",
        "lambu = wishart.rvs(size = 1,df = vu0, scale = wu0)\n",
        "muu = np.random.multivariate_normal(muu0, np.linalg.pinv(beta0 * lambu))\n",
        "muv = np.random.multivariate_normal(muv0, np.linalg.pinv(beta0 * lambv))\n",
        "U = np.zeros((N,D))  \n",
        "V = np.zeros((M,D))\n",
        "U = np.random.multivariate_normal(muu,lambu,N)\n",
        "V = np.random.multivariate_normal(muv,lambv,M)\n",
        "sigmaUN = np.zeros((N, D, D))\n",
        "sigmaVM = np.zeros((M, D, D))\n",
        "muun = np.zeros((N, D))\n",
        "muvm = np.zeros((M, D))\n",
        "vvt = np.zeros((M , D, D))\n",
        "uut = np.zeros((N , D, D))\n",
        "for i in range(M):\n",
        "  vvt[i] = V[i, :].reshape(D, 1) @ V[i ,:].reshape(1,D)\n",
        "for i in range(N):\n",
        "  uut[i] = U[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)  \n",
        " "
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aJWy4IKL6dzl",
        "colab_type": "code",
        "outputId": "5ee05598-5e20-4ee3-ad28-473dcf290263",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1000
        }
      },
      "source": [
        "num = 100\n",
        "U_spl = np.zeros((num ,N , D))\n",
        "V_spl = np.zeros((num ,M ,D))\n",
        "for t in range(num):\n",
        "  \n",
        "  #Calculating joint posterior for lambdaU and muU\n",
        "  \n",
        "  ubar = np.sum(U ,axis = 0) / N\n",
        "  #print(ubar,\"bar\")\n",
        "  subar = np.zeros((D ,D))\n",
        "  for i in range(N):\n",
        "    subar += U[i ,:].reshape(D ,1) @ U[i ,:].reshape(1 ,D)\n",
        "  subar /= N\n",
        "  #print(subar,\"subar\")\n",
        "  betastaru = beta0 + N\n",
        "  muustar = (beta0 * mu0 + N * ubar) / betastaru  \n",
        "  #print(muustar,\"!!!!\")\n",
        "  vustar = v0 + N\n",
        "  wustr0 = np.linalg.pinv(np.linalg.pinv(w0) + N * subar + beta0 * N * (mu0 - ubar).reshape(D ,1) @ (mu0 - ubar).reshape(1 ,D) / betastaru)\n",
        "  \n",
        " #Calculating joint posterior for lambdaV and muV\n",
        "  \n",
        "  vbar = np.sum(V ,axis = 0) / M\n",
        "  svbar = np.zeros((D , D))\n",
        "  for j in range(M):\n",
        "    svbar += V[i ,:].reshape(D ,1) @ V[i ,:].reshape(1 ,D)\n",
        "  svbar /= M\n",
        "  muvstar = (beta0 * mu0 + M * vbar) / (beta0 + M) #muv0 pass\n",
        "  #print(muvstar,\"1234\")\n",
        "  betastrv = beta0 + M\n",
        "  vvstar = v0 + M\n",
        "  wvstr0 = np.linalg.pinv(np.linalg.pinv(w0) + M * svbar + beta0 * M * (mu0 - vbar).reshape(D ,1) @ (mu0 - vbar).reshape(1 ,D) / betastrv)\n",
        " \n",
        "  # Sampling lambaU and muU\n",
        "  \n",
        "  lambu = wishart.rvs(size = 1,df = vustar, scale = wustr0)\n",
        "  muu = np.random.multivariate_normal(muustar, np.linalg.pinv(betastaru * lambu))\n",
        "  #print(lambu,\"lambu\")\n",
        "  #print(muu,\"mu\")\n",
        "   \n",
        "  # Sampling lambdaV and muV\n",
        "  \n",
        "  lambv = wishart.rvs(size = 1,df = vvstar, scale = wvstr0)\n",
        "  muv = np.random.multivariate_normal(muvstar, np.linalg.pinv(betastrv * lambv))\n",
        "  #print(muv,\"muv\")\n",
        "  \n",
        "  #Calculating Conditional Posterior for U\n",
        "  \n",
        "  for i in range(N):\n",
        "    s1 = np.zeros((D, D))\n",
        "    s2 = np.zeros(D)\n",
        "    for j in range(M):\n",
        "      s1 += I[i][j] * vvt[j]\n",
        "      s2 += R[i][j] * V[j ,:]\n",
        "    s2 *= alpha\n",
        "    sigmaUN[i] = np.linalg.pinv(lambu + alpha * s1)\n",
        "    muun[i] = sigmaUN[i] @ ( s2 + lambu @ muu)\n",
        "    \n",
        "  #Sampling U\n",
        "  \n",
        "  for i in range(N):\n",
        "    U[i] = np.random.multivariate_normal(muun[i], sigmaUN[i])\n",
        "  U_spl[t] = U\n",
        "  #print(U_spl[t])\n",
        "  for i in range(N):\n",
        "    uut[i] = U[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)\n",
        "    \n",
        "  #Calculating Conditional Posterior for V\n",
        "  \n",
        "  for j in range(M):#\n",
        "    A = np.zeros((D ,D))\n",
        "    B = np.zeros(D)\n",
        "    for i in range(N):\n",
        "      A += I[i][j] * uut[i]\n",
        "      B += R[i][j] * U[i ,:]\n",
        "    B *= alpha\n",
        "    sigmaVM[j] = np.linalg.pinv(lambv + alpha * A)\n",
        "    muvm[j] = sigmaVM[j] @ (B + lambv @ muv)\n",
        "    \n",
        "  #Sampling V\n",
        "  \n",
        "  for j in range(M):\n",
        "    V[j] = np.random.multivariate_normal(muvm[j] ,sigmaVM[j])\n",
        "  V_spl[t] = V\n",
        "  #print(V_spl[t])\n",
        "  for i in range(M):\n",
        "    vvt[i] = V[i, :].reshape(D, 1) @ V[i ,:].reshape(1, D)\n",
        "    \n",
        "  print (t + 1)\n",
        "  predictR = np.zeros((N,M))\n",
        "  predictR += U_spl[t]@V_spl[t].T\n",
        "  #print(predictR)\n",
        "  RMSE = (np.sum(((ITest*predictR) - RTest) ** 2) / np.sum(ITest)) ** 0.5\n",
        "  print(RMSE)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "1\n",
            "5.920613308822863\n",
            "2\n",
            "1.4707984611868525\n",
            "3\n",
            "1.0231974995193707\n",
            "4\n",
            "1.0150364654538937\n",
            "5\n",
            "1.0156753062154362\n",
            "6\n",
            "1.0288489979147588\n",
            "7\n",
            "1.0217198500795015\n",
            "8\n",
            "1.0145596372864243\n",
            "9\n",
            "1.0278115957983862\n",
            "10\n",
            "1.0293119925407916\n",
            "11\n",
            "1.0416835858029616\n",
            "12\n",
            "1.033870651584277\n",
            "13\n",
            "1.0373943327482311\n",
            "14\n",
            "1.0360493648450388\n",
            "15\n",
            "1.0447931088026012\n",
            "16\n",
            "1.0506483427840043\n",
            "17\n",
            "1.045830370104952\n",
            "18\n",
            "1.0429110683886613\n",
            "19\n",
            "1.0418252931768759\n",
            "20\n",
            "1.0509166957937979\n",
            "21\n",
            "1.054067451159698\n",
            "22\n",
            "1.045856464566836\n",
            "23\n",
            "1.0488898986638941\n",
            "24\n",
            "1.0532168422408785\n",
            "25\n",
            "1.0409634191145367\n",
            "26\n",
            "1.0504744446020953\n",
            "27\n",
            "1.0520367953513543\n",
            "28\n",
            "1.053202278311832\n",
            "29\n",
            "1.0581523503623895\n",
            "30\n",
            "1.0495016644441029\n",
            "31\n",
            "1.052708003379424\n",
            "32\n",
            "1.0594354937194164\n",
            "33\n",
            "1.0570331339844274\n",
            "34\n",
            "1.0519092389614924\n",
            "35\n",
            "1.0516170314466144\n",
            "36\n",
            "1.0585206105888776\n",
            "37\n",
            "1.0433409934097746\n",
            "38\n",
            "1.0619878894441754\n",
            "39\n",
            "1.0533186026335335\n",
            "40\n",
            "1.0571418975847358\n",
            "41\n",
            "1.0472438764331586\n",
            "42\n",
            "1.0487649778831536\n",
            "43\n",
            "1.0526607120357576\n",
            "44\n",
            "1.057935897662006\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AR-YIe8aS9bm",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "_71oFsPSr8Sr",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "predictR = np.zeros((N,M))\n",
        "for i in range(num):\n",
        "    predictR += U_spl[i]@V_spl[i].T\n",
        "predictR/=num\n",
        "predCheck = ITest * predictR\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "MZLFB53pA1jW",
        "colab_type": "code",
        "outputId": "adabcb42-f488-4275-98d5-2cee2d16171b",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 33
        }
      },
      "source": [
        "RMSE=0\n",
        "RMSE = (np.sum((predCheck - RTest) ** 2) / np.sum(ITest)) ** 0.5\n",
        "RMSE"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "0.9197090006938149"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 21
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "oBKtp4Moqs6R",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "co_Lcb4wRMEm",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def sample_mu(size,uv,lamb,mu):\n",
        "    bar = np.sum(uv ,axis = 0) / size\n",
        "    #print(bar,\"bar\")\n",
        "    sbar = np.zeros((D ,D))\n",
        "    for i in range(size):\n",
        "        sbar += uv[i ,:].reshape(D ,1) @ uv[i ,:].reshape(1 ,D) \n",
        "    sbar /= size\n",
        "    #print(subar,\"subar\")\n",
        "    mustar = (beta0 * muu0 + size * bar) / (beta0 + size)\n",
        "    #print(mustar,\"muuuu\")\n",
        "    betastar = beta0 + size\n",
        "    vstar = vu0 + size\n",
        "    wstr0 = np.linalg.pinv(np.linalg.pinv(wu0) + N * sbar + beta0*size*(muu0 - bar).reshape(D ,1) @ (muu0 - bar).reshape(1 ,D) / betastar)\n",
        "\n",
        "    lamb = wishart.rvs(size = 1,df = vstar, scale = wstr0)\n",
        "    #print(lamb,\"lambu\")\n",
        "    mu = np.random.multivariate_normal(mustar, np.linalg.pinv(betastar * lamb))\n",
        "    print(mu,\"mu\")\n",
        "\n",
        "    return lamb,mu\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "DRhlqfGYAFCU",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def sample_u(lamb,mu,uv,uv2,size,size2):\n",
        "    sigma = np.zeros((size, D, D))\n",
        "    mun = np.zeros((size, D))\n",
        "    for i in range(size): \n",
        "      s1 = np.zeros((D,D))\n",
        "      s2 = np.zeros(D)\n",
        "      for j in range(size2):\n",
        "        s1 += I[i][j] * uv2[i, :].reshape(D, 1) @ uv2[i ,:].reshape(1,D)#vut[j]\n",
        "        s2 += R[i][j] * uv2[j ,:]#uv2[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)\n",
        "#uv2[j ,:]\n",
        "      \n",
        "      s2 *= alpha\n",
        "      sigma[i] = np.linalg.pinv(lamb + alpha * s1)\n",
        "      mun[i] = sigma[i] @ ( s2 + lamb @ mu)\n",
        "    \n",
        "    #Sampling U\n",
        "  \n",
        "    for i in range(size):\n",
        "      uv[i] = np.random.multivariate_normal(mun[i], sigma[i])\n",
        "    spl = uv\n",
        "    print(spl[t])\n",
        "    # for i in range(size):\n",
        "    #   uut[i] = U[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)\n",
        "    return spl\n",
        "def sample_v(lamb,mu,uv,uv2,size,size2):\n",
        "    sigma = np.zeros((size, D, D))\n",
        "    mun = np.zeros((size, D))\n",
        "    for j in range(size): \n",
        "      s1 = np.zeros((D,D))\n",
        "      s2 = np.zeros(D)\n",
        "      for i in range(size2):\n",
        "        s1 += I[i][j] * uv2[i, :].reshape(D, 1) @ uv2[i ,:].reshape(1,D)#vut[j]\n",
        "        s2 += R[i][j] * uv2[i ,:]#uv2[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)\n",
        "#uv2[j ,:]\n",
        "      \n",
        "      s2 *= alpha\n",
        "      sigma[j] = np.linalg.pinv(lamb + alpha * s1)\n",
        "      mun[j] = sigma[j] @ ( s2 + lamb @ mu)\n",
        "    \n",
        "    #Sampling U\n",
        "  \n",
        "    for i in range(size):\n",
        "      uv[i] = np.random.multivariate_normal(mun[i], sigma[i])\n",
        "    spl = uv\n",
        "    print(spl[t])\n",
        "    # for i in range(size):\n",
        "    #   uut[i] = U[i, :].reshape(D, 1) @ U[i ,:].reshape(1,D)\n",
        "    return spl\n"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ZFbd9LerkdsF",
        "colab_type": "code",
        "outputId": "74355734-cf40-4301-88f8-9241cb217596",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 375
        }
      },
      "source": [
        "T = 100\n",
        "U_spl = np.zeros((T ,N , D))\n",
        "V_spl = np.zeros((T ,M ,D))\n",
        "for t in range(T):\n",
        "\n",
        "  #sample log and mu\n",
        "  lambu,muu = sample_mu(N,U,lambu,muu)\n",
        "  #print(lambu)\n",
        "  #sample log and mu\n",
        "  lambv,muv = sample_mu(M,V,lambv,muv)\n",
        "\n",
        "  #sample U\n",
        "  U_spl[t] = sample_u(lambu,muu,U,V,N,M)\n",
        "  V_spl[t] = sample_v(lambv,muv,V,U,M,N)\n",
        "  print (t + 1)\n",
        "  predictR = np.zeros((N,M))\n",
        "  predictR += U_spl[t]@V_spl[t].T\n",
        "  print(predictR)\n",
        "  RMSE = (np.sum(((ITest*predictR) - RTest) ** 2) / np.sum(ITest)) ** 0.5\n",
        "  print(RMSE)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "[-1.4713551  -1.4561111   0.53010859  0.08178807  0.02132132  0.806325\n",
            " -0.56314941  1.67429604 -0.20121779 -0.77969724 -0.43481091 -0.31661084\n",
            "  1.68595379 -0.89252587 -0.98471954  0.17153333  0.9594806   0.57905774\n",
            "  1.43526355 -0.24766978 -0.17960006 -1.41602098  0.90822037  2.85659318\n",
            "  1.64856462 -1.40181524  1.9469802  -1.80814997 -1.79881124  3.56562302\n",
            "  0.10019048 -1.24979963  0.21468621  0.71219619 -0.72706254  0.88993962\n",
            "  0.39381166  1.07104053 -0.80185257  0.92278961  0.10321413  3.0389393\n",
            " -0.11418408  0.79478127  1.8502203  -0.1742113  -2.22596288  1.07462453\n",
            "  0.05319862  0.75414058] mu\n",
            "[-5.59804548e-01  2.34416497e+00 -5.61930136e+00 -5.19482604e+00\n",
            "  9.39808913e+00  5.37966750e+00 -8.40259866e-01  1.03218114e+00\n",
            "  6.45461059e+00  5.10153946e+00  2.07788653e+00  1.42499687e+00\n",
            "  3.76645558e+00  6.35869728e-01 -3.49744039e+00 -2.03550166e+00\n",
            "  3.07098222e-02  2.58164927e+00 -3.01824133e+00  1.98360676e+00\n",
            "  3.65087784e+00 -2.48085448e-01  2.00767013e+00  7.49137761e+00\n",
            "  4.00835332e+00 -9.52742131e-01  7.47608539e-01  8.82691982e-01\n",
            "  6.53602941e+00  2.33267344e+00  2.75569330e+00  8.95993308e+00\n",
            "  1.10477078e+00  1.41987334e-01 -9.53822710e-01  6.20484908e+00\n",
            "  5.41915121e+00 -1.98275876e+00  7.21230334e+00 -1.51715611e+00\n",
            " -2.80798753e+00  5.14153907e-01  2.35789694e-01  3.94260724e-03\n",
            " -2.13690689e+00 -5.39816223e-01  4.03833649e+00  8.02637821e-01\n",
            " -5.15908905e+00 -2.67495541e-01] mu\n"
          ],
          "name": "stdout"
        }
      ]
    }
  ]
}