pgmpy_intro.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "pgmpy.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyNJ0KNWO5XTYNQdk5r8B4Jw",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/kpapdac/graphicalModels/blob/main/pgmpy_intro.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Wb_hz1foxg98",
        "outputId": "124b8dac-c651-428a-b268-3dcc9b4e7872"
      },
      "source": [
        "%pip install pgmpy"
      ],
      "execution_count": 1,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Collecting pgmpy\n",
            "  Downloading pgmpy-0.1.16-py3-none-any.whl (1.9 MB)\n",
            "\u001b[K     |████████████████████████████████| 1.9 MB 6.5 MB/s \n",
            "\u001b[?25hRequirement already satisfied: tqdm in /usr/local/lib/python3.7/dist-packages (from pgmpy) (4.62.3)\n",
            "Requirement already satisfied: numpy in /usr/local/lib/python3.7/dist-packages (from pgmpy) (1.19.5)\n",
            "Requirement already satisfied: pandas in /usr/local/lib/python3.7/dist-packages (from pgmpy) (1.1.5)\n",
            "Requirement already satisfied: scipy in /usr/local/lib/python3.7/dist-packages (from pgmpy) (1.4.1)\n",
            "Requirement already satisfied: networkx in /usr/local/lib/python3.7/dist-packages (from pgmpy) (2.6.3)\n",
            "Requirement already satisfied: torch in /usr/local/lib/python3.7/dist-packages (from pgmpy) (1.9.0+cu111)\n",
            "Requirement already satisfied: scikit-learn in /usr/local/lib/python3.7/dist-packages (from pgmpy) (0.22.2.post1)\n",
            "Requirement already satisfied: joblib in /usr/local/lib/python3.7/dist-packages (from pgmpy) (1.0.1)\n",
            "Requirement already satisfied: statsmodels in /usr/local/lib/python3.7/dist-packages (from pgmpy) (0.10.2)\n",
            "Requirement already satisfied: pyparsing in /usr/local/lib/python3.7/dist-packages (from pgmpy) (2.4.7)\n",
            "Requirement already satisfied: pytz>=2017.2 in /usr/local/lib/python3.7/dist-packages (from pandas->pgmpy) (2018.9)\n",
            "Requirement already satisfied: python-dateutil>=2.7.3 in /usr/local/lib/python3.7/dist-packages (from pandas->pgmpy) (2.8.2)\n",
            "Requirement already satisfied: six>=1.5 in /usr/local/lib/python3.7/dist-packages (from python-dateutil>=2.7.3->pandas->pgmpy) (1.15.0)\n",
            "Requirement already satisfied: patsy>=0.4.0 in /usr/local/lib/python3.7/dist-packages (from statsmodels->pgmpy) (0.5.2)\n",
            "Requirement already satisfied: typing-extensions in /usr/local/lib/python3.7/dist-packages (from torch->pgmpy) (3.7.4.3)\n",
            "Installing collected packages: pgmpy\n",
            "Successfully installed pgmpy-0.1.16\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eXxnoDvc7f1x"
      },
      "source": [
        "### CPD for Discrete Variables"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fu74IZFZ1ec7"
      },
      "source": [
        "Create network structure"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "FQrRHdYKynOo"
      },
      "source": [
        "from pgmpy.models import BayesianNetwork\n",
        "from pgmpy.factors.discrete.CPD import TabularCPD\n",
        "import pandas as pd\n",
        "import numpy as np\n",
        "model = BayesianNetwork([('A', 'B'), ('B', 'C'), ('A', 'D'), ('D', 'C')])"
      ],
      "execution_count": 6,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "unyH4A1i2Fcn"
      },
      "source": [
        "Create CPDs from dataframe"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "geGKWjAh163N",
        "outputId": "cf1a4788-898b-4e0c-872c-c9e88406f9de"
      },
      "source": [
        "values = pd.DataFrame(np.random.randint(low=0, high=2, size=(1000, 4)), columns=['A', 'B', 'C', 'D'])\n",
        "model.fit(values)\n",
        "model.get_cpds()"
      ],
      "execution_count": 7,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[<TabularCPD representing P(A:2) at 0x7f00ed7b3410>,\n",
              " <TabularCPD representing P(B:2 | A:2) at 0x7f00ef08ded0>,\n",
              " <TabularCPD representing P(C:2 | B:2, D:2) at 0x7f00ef028f10>,\n",
              " <TabularCPD representing P(D:2 | A:2) at 0x7f00ef08df50>]"
            ]
          },
          "metadata": {},
          "execution_count": 7
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eV3nJbYU2OB0"
      },
      "source": [
        "Check model validity"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "0vIuscsb2Mhb",
        "outputId": "fe57c192-4093-452b-a992-41b2de2cb109"
      },
      "source": [
        "model.check_model()"
      ],
      "execution_count": 8,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "True"
            ]
          },
          "metadata": {},
          "execution_count": 8
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Du8VnDzn2Y11"
      },
      "source": [
        "Manually assigning CPDs"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "_u7F50VV2V0n",
        "outputId": "58f09dfd-5957-4edb-e6c1-5deaa175cc9f"
      },
      "source": [
        "model = BayesianNetwork([('D', 'G'), ('I', 'G'), ('G', 'L'), ('I', 'S')])\n",
        "# Defining individual CPDs.\n",
        "cpd_d = TabularCPD(variable='D', variable_card=2, values=[[0.6], [0.4]])\n",
        "cpd_i = TabularCPD(variable='I', variable_card=2, values=[[0.7], [0.3]])\n",
        "\n",
        "# The representation of CPD in pgmpy is a bit different than the CPD shown in the above picture. In pgmpy the colums\n",
        "# are the evidences and rows are the states of the variable. So the grade CPD is represented like this:\n",
        "#\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | diff    | intel_0 | intel_0 | intel_1 | intel_1 |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | intel   | diff_0  | diff_1  | diff_0  | diff_1  |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_0 | 0.3     | 0.05    | 0.9     | 0.5     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_1 | 0.4     | 0.25    | 0.08    | 0.3     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_2 | 0.3     | 0.7     | 0.02    | 0.2     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "\n",
        "cpd_g = TabularCPD(variable='G', variable_card=3,\n",
        "                   values=[[0.3, 0.05, 0.9,  0.5],\n",
        "                           [0.4, 0.25, 0.08, 0.3],\n",
        "                           [0.3, 0.7,  0.02, 0.2]],\n",
        "                  evidence=['I', 'D'],\n",
        "                  evidence_card=[2, 2])\n",
        "\n",
        "cpd_l = TabularCPD(variable='L', variable_card=2,\n",
        "                   values=[[0.1, 0.4, 0.99],\n",
        "                           [0.9, 0.6, 0.01]],\n",
        "                   evidence=['G'],\n",
        "                   evidence_card=[3])\n",
        "\n",
        "cpd_s = TabularCPD(variable='S', variable_card=2,\n",
        "                   values=[[0.95, 0.2],\n",
        "                           [0.05, 0.8]],\n",
        "                   evidence=['I'],\n",
        "                   evidence_card=[2])\n",
        "\n",
        "# Associating the CPDs with the network\n",
        "model.add_cpds(cpd_d, cpd_i, cpd_g, cpd_l, cpd_s)\n",
        "\n",
        "# check_model checks for the network structure and CPDs and verifies that the CPDs are correctly\n",
        "# defined and sum to 1.\n",
        "model.check_model()"
      ],
      "execution_count": 9,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "True"
            ]
          },
          "metadata": {},
          "execution_count": 9
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "04dZwtC07k1E"
      },
      "source": [
        "### CPD for Continuous Variables"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TZkybh2L0on6"
      },
      "source": [
        "Dirichlet distribution of two continuous variables x and y"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "SBi7cPk-218I",
        "outputId": "6d3c48e7-e03c-4f13-c950-c83097d0d654"
      },
      "source": [
        "import numpy as np\n",
        "from scipy.special import beta\n",
        "\n",
        "# Two variable dirichlet ditribution with alpha = (1,2)\n",
        "def dirichlet_pdf(x, y):\n",
        "     return (np.power(x, 1)*np.power(y, 2))/beta(1,2)\n",
        "\n",
        "from pgmpy.factors.continuous import ContinuousFactor\n",
        "dirichlet_factor = ContinuousFactor(['x', 'y'], dirichlet_pdf)\n",
        "dirichlet_factor.scope(), dirichlet_factor.assignment(5,6)"
      ],
      "execution_count": 42,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "(['x', 'y'], 360.0)"
            ]
          },
          "metadata": {},
          "execution_count": 42
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 283
        },
        "id": "kJqw0swn22u5",
        "outputId": "a6f11d89-4550-41f6-cf70-454eec911383"
      },
      "source": [
        "import matplotlib.pyplot as plt\n",
        "plt.plot(dirichlet_pdf(np.linspace(0,1,100),1))"
      ],
      "execution_count": 38,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[<matplotlib.lines.Line2D at 0x7f00e7703710>]"
            ]
          },
          "metadata": {},
          "execution_count": 38
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deXgV1f3H8feX7IFAIESWLARkk1UgbGqtbUVxxaUWVBQoiLa1Wlt/rXaxrdrN2tpardYFEREQrW0RUeqCiguQALKFLUQgYcvCGkLWe35/5OqTIpAAN5ncez+v58nDnblD5jvMzYfJmTnnmHMOEREJfi28LkBERAJDgS4iEiIU6CIiIUKBLiISIhToIiIhItKrHbdv395lZGR4tXsRkaC0fPnyYudc8rHe8yzQMzIyyM7O9mr3IiJBycy2He89NbmIiIQIBbqISIhQoIuIhAgFuohIiFCgi4iEiHoD3cymmVmhma09zvtmZo+aWa6ZrTazwYEvU0RE6tOQK/TpwOgTvH8J0MP/NRV44vTLEhGRk1VvoDvnPgD2nmCTMcAMV2sJkGhmnQJVoIhIqPD5HL95PYf8vWWN8v0D0YaeAuTXWS7wr/sSM5tqZtlmll1UVBSAXYuIBI+/vZvL04s/48Pc4kb5/k16U9Q595RzLtM5l5mcfMyeqyIiIen9TUX85Z1NXDM4hXFD0xplH4EI9B1A3epS/etERATYsf8Id85ZSa8OCfzmqv6YWaPsJxCBPg+42f+0ywjggHNuVwC+r4hI0KuoruG7L66gpsbxxPghxEVHNNq+6h2cy8xmAxcA7c2sAPglEAXgnHsSWABcCuQCZcCkxipWRCTY3P9aDqvy9/Pk+CF0bd+yUfdVb6A7566v530HfC9gFYmIhIiXs/N5cel2bvvqmYzu17HR96eeoiIijWDtjgP87N9rObd7Endf1LNJ9qlAFxEJsP1lldw2czlJLaN5dNwgIiOaJmo9m+BCRCQU1fgcd8z5lMKDFbx06wiSWsU02b4V6CIiAfTIW5v4YFMRv7umP4PS2zbpvtXkIiISIAvX7eaxRbmMG5rG9cPSm3z/CnQRkQDILSzlR3NXMTC1Db+6sq8nNSjQRURO06HyKqa+kE1MZAueGD+E2KjG6zx0ImpDFxE5DT6f44dzV7GtpIwXpwync2KcZ7XoCl1E5DQ8tiiXt3L28PPLzmJEtyRPa1Ggi4iconfW7+GRtzdxzaAUJp6T4XU5CnQRkVORW1jKD+Z8St/OrfntNY03guLJUKCLiJykg/6boNGRLfjHTZme3QQ9mm6KioicBJ/P8YM5n7LdfxM0xcOboEfTFbqIyEn481ubeHdDIb+8og/DPb4JejQFuohIA81fvfOLnqDjR3TxupwvUaCLiDTA2h0HuPvlVWR2acv9Y/o1i5ugR1Ogi4jUo7i0gqkzsmkbH80T44cQHdk8o1M3RUVETqCy2sd3Zi6n5HAlr9x2DskJTTcc7slSoIuIHIdzjl/8ey1ZW/fxt+sH0T+1jdclnVDz/L1BRKQZmP7xVl7Kzuf2r3XnioGdvS6nXgp0EZFjWLy5iAfm53BRnw78cFTTzAl6uhToIiJH2VJUyndfXEHPDgk8MvZsWrRofk+0HIsCXUSkjv1llUx5PpvoiBY8fXMmLWOC51Zj8FQqItLIqmp8fG/WCnbsO8KLtwwnrV281yWdFAW6iAi1T7T8+rV1fJRbwsPXDWRoRjuvSzppanIREQGe/3grM5ds59bzu/HNIalel3NKFOgiEvYWbSzk/vk5jOrTgZ+M7u11OadMgS4iYW3j7kN8f9ZKendszV+C6ImWY1Ggi0jYKjpUweTns4iPjuDZicH1RMuxBHf1IiKnqLyqhqkvZFNcWsHcW0fSqU3zmajiVCnQRSTs+HyOu19excrt+3ly/GAGpCZ6XVJANKjJxcxGm9lGM8s1s3uO8X66mS0ys5VmttrMLg18qSIigfHI25uYv3oX91zSm9H9OnldTsDUG+hmFgE8DlwC9AGuN7M+R232c2Cuc24QMA74e6ALFREJhJez8/nbu7mMzUzj1vO7eV1OQDXkCn0YkOucy3POVQJzgDFHbeOA1v7XbYCdgStRRCQwPs4t5t5X13Bu9yQevLp5zjp0OhoS6ClAfp3lAv+6un4FjDezAmAB8P1jfSMzm2pm2WaWXVRUdArlioicmtzCQ9w6czld27fk7zcOISoi9B7yC9QRXQ9Md86lApcCL5jZl763c+4p51ymcy4zOTk5QLsWETmxokMVTHwui5jICKZNHEqbuCivS2oUDQn0HUBaneVU/7q6JgNzAZxznwCxQPtAFCgicjrKKquZ8nwWJaWVPDshM+gG3DoZDQn0LKCHmXU1s2hqb3rOO2qb7cA3AMzsLGoDXW0qIuKpGp/jjtkrWbPjAI9eP4iBaaHxeOLx1Bvozrlq4HZgIbCe2qdZ1pnZ/WZ2pX+zHwG3mNkqYDYw0TnnGqtoEZH6OOe4/7V1vL2+kF9e0ZdRfTp4XVKja1DHIufcAmpvdtZdd1+d1znAuYEtTUTk1D2z+DOe/2QbU87ryoRzMrwup0mE3m1eEQl781bt5DcL1nPZgE789NKzvC6nySjQRSSkLMkr4e65qxiW0Y4/XTcwqEdPPFkKdBEJGZv3HGLqjGzS2sXx1M1DiI2K8LqkJqVAF5GQsPtAOROmLSMmKoLpk4aRGB/tdUlNToEuIkHvYHkVE59bxsHyaqZPGhrSz5qfiAJdRIJaRXUNt85YTm5hKU+OH0Lfzm28LskzGg9dRIKWz+f44dxVfJJXwiNjB3Jej/DuoK4rdBEJSs45fv3aOl5fvYufXtqbqwelel2S5xToIhKU/v7eli86Dk09/0yvy2kWFOgiEnTmZuXzx4UbuerszmHVcag+CnQRCSoL1+3mnldXc37PZB76Znh1HKqPAl1EgsbSvBK+P3sl/VMTeeLGwURHKsLq0r+GiASFnJ0HmTIjm7S2cTw3cSgtY/SQ3tEU6CLS7G0tPszN05bRKiaSGZOH065l+PUCbQgFuog0a3sOljP+2aXU+Hy8MHkYKYlxXpfUbCnQRaTZOlBWxc3PLmPf4UqmTxpG9zMSvC6pWVMjlIg0S4crqpk4fRmfFR/muUlDQ376uEBQoItIs1NRXcPUF7JZXXCAx28YzLndw7tLf0OpyUVEmpXqGh93zF7JR7klPHTtAEb36+h1SUFDgS4izYbP5/jxP1ezcN0efnlFH64dovFZToYCXUSaBeccv3ptHa+u2MFdF/Zk0rldvS4p6CjQRaRZeGjhRmZ8so2p53fjjm9097qcoKRAFxHPPb4olyfe28KNw9O595LemGl8llOhQBcRTz374Wf8ceFGrh6UwgNj+inMT4MCXUQ8M2vpdh6Yn8Ml/Tryx28O0MiJp0mBLiKeeHVFAT/79xq+1iuZv44bRGSE4uh06V9QRJrc/NU7ufvlVYzslsQT44doGNwA0b+iiDSphet2c+ecT8ns0o5nJmQSGxXhdUkhQ4EuIk1m0YZCbp+1ggGpbZg2aSjx0Rp9JJAU6CLSJN7fVMStM5fTq2MC0ycNo5UmqAi4BgW6mY02s41mlmtm9xxnm2+ZWY6ZrTOzWYEtU0SC2Yebi5k6I5vuya2YOXk4beKivC4pJNX7X6SZRQCPA6OAAiDLzOY553LqbNMDuBc41zm3z8zOaKyCRSS4fLylmCkzsujaviUzpwwnMV6zDTWWhlyhDwNynXN5zrlKYA4w5qhtbgEed87tA3DOFQa2TBEJRkvySpg8PZv0dvG8OEVTxzW2hgR6CpBfZ7nAv66unkBPM/vIzJaY2ehjfSMzm2pm2WaWXVRUdGoVi0hQWJJXwqTnskhtG8eLU0aQ1CrG65JCXqBuikYCPYALgOuBp83sS9OLOOeecs5lOucyk5OTA7RrEWlulvrDPKVtHLNuGUFygsK8KTQk0HcAaXWWU/3r6ioA5jnnqpxznwGbqA14EQkzS/JKmPhcFp0TY5l1y3CFeRNqSKBnAT3MrKuZRQPjgHlHbfNvaq/OMbP21DbB5AWwThEJAh9vKf6imWXO1JGckRDrdUlhpd5Ad85VA7cDC4H1wFzn3Dozu9/MrvRvthAoMbMcYBHwf865ksYqWkSan49yi/n29CzS2qmZxSvmnPNkx5mZmS47O9uTfYtIYL2/qYipM7LJSGrJi7cMp71ugDYaM1vunMs81nvqqiUip+Wd9Xv4zswVnHlGK2ZOHqanWTykQBeRU7Zw3W5un7WC3h1b88LkYeo05DEFuoicktdW7eSulz6lX0obnv/2MHXnbwY0OJeInLRXlhdw55yVDEpP5IXJCvPmQlfoInJSZi3dzk//tYZzuyfx9M2ZGgK3GdGZEJEGe2ZxHg++vp6v9UrmifFDNDlFM6NAF5F6Oef427u5/PmtTVzSryN/HTdI08Y1Qwp0ETkh5xy/f3MD/3g/j2sGp/DQtQM0oXMzpUAXkeOq8Tl+8Z+1zFq6nfEj0rn/yn60aGFelyXHoUAXkWOqqvFx98ur+M+nO/nOBWfy44t7YaYwb84U6CLyJeVVNdw+awVvry/kx6N78d0LuntdkjSAAl1E/sfB8iqmPJ9N1ta9PHBVP24a0cXrkqSBFOgi8oXi0gomTFvGxt2H+MvYsxlz9tGTk0lzpkAXEQAK9pVx87PL2HngCM9MyOSCXprrPdgo0EWEDbsPMmHaMo5U1vDilOEM6dLO65LkFCjQRcJc1ta9TJ6eRVx0BC/fdg69OiZ4XZKcIgW6SBh7O2cP35u1gpTEOGZMHkZq23ivS5LToEAXCVNzltUOstUvpQ3PTRyqiSlCgAJdJMzUHZflqz2T+fuNg2kZoygIBTqLImGkusbHffPWMWvpdq4ZnMIfrh1AlMZlCRkKdJEwUVZZzR2zV/L2+kJ15Q9RCnSRMFBSWsG3n89mTcF+HhjTl5tGZnhdkjQCBbpIiMsrKmXS9Cx2HyjnyfFDuKhvR69LkkaiQBcJYdlb9zJlRjYtzJg9dQSD09t6XZI0IgW6SIh6ffUu7pr7KSmJcUyfNJQuSS29LkkamQJdJMQ453ji/S089OZGMru05embM2nbMtrrsqQJKNBFQkhVjY+f/2stL2Xnc+XAzjz0zQGayDmMKNBFQsSBsiq+O2s5H+WWcMfXu3PXqJ56LDHMKNBFQsC2ksNMmp5F/t4yHr5uIN8ckup1SeIBBbpIkFuaV8JtM5fjgJmThzO8W5LXJYlHFOgiQWxuVj4/+/ca0trFM23CUDLa60mWcNagQRzMbLSZbTSzXDO75wTbXWtmzswyA1eiiBytxud4cH4OP/7nakZ0S+Jf3z1XYS71X6GbWQTwODAKKACyzGyecy7nqO0SgDuBpY1RqIjUOlhexZ2zV7JoYxETRnbhF5f3IVIDbAkNu0IfBuQ65/Kcc5XAHGDMMbZ7APgDUB7A+kSkjs+KD3P14x+xeHMxD1zVj1+P6acwly805JOQAuTXWS7wr/uCmQ0G0pxzr5/oG5nZVDPLNrPsoqKiky5WJJwt3lzEmMc+ZO/hSl6YPJybRnTxuiRpZk77v3YzawH8GfhRfds6555yzmU65zKTk5NPd9ciYcE5x9Mf5DFh2jI6tYlj3u3nMfJMPckiX9aQp1x2AGl1llP96z6XAPQD3vN3YugIzDOzK51z2YEqVCQclVfVcO+ra/jXyh2M7tuRP31roGYXkuNqyCcjC+hhZl2pDfJxwA2fv+mcOwC0/3zZzN4D7laYi5yeHfuPcNsLy1m78wB3X9ST732tu3p+ygnVG+jOuWozux1YCEQA05xz68zsfiDbOTevsYsUCTcf5xZz++yVVFX7ePqmTC7s08HrkiQINOh3N+fcAmDBUevuO862F5x+WSLhyTnHM4s/43dvrOfM5Fb846YhdEtu5XVZEiTUGCfSTJRWVPOTV1bz+ppdXNKvI3+8biCt1F4uJ0GfFpFmILewlNtmLievqJR7L+nN1PO7qb1cTpoCXcRj81fv5CevrCY2KoKZU4Zzzpnt6/9LIsegQBfxSGW1j98uWM/0j7cyOD2Rx28cTKc2cV6XJUFMgS7igZ37j/C9WStYuX0/3z63K/de2psodeGX06RAF2liizYUctfcT6mucTx+w2AuG9DJ65IkRCjQRZpIdY2Ph/+7iSff30KfTq15/MbBdNWQtxJACnSRJrBj/xHumL2S5dv2ccPwdO67vI8mb5aAU6CLNLK3cvZw98urqK7x8ej1g7hyYGevS5IQpUAXaSQV1TX8bsEGpn+8lb6dW/PYDWpikcalQBdpBLmFpdwxeyU5uw4y8ZwM7r20NzGRamKRxqVAFwkg5xwvZeXz69dyiI1qwbMTMvnGWRpYS5qGAl0kQPaXVXLvq2t4Y+1uRnZL4pGxZ9OxTazXZUkYUaCLBMDHW4r54UurKDlcwb2X9OaWr3SjRQuNxSJNS4Euchoqqmv403838fTiPLomteTpm8+lf2obr8uSMKVAFzlFG3cf4s45K9mw+xDjR6Tz00vPIj5aP1LiHX36RE5Sjc/x7Id5PLxwE63jIpk2MZOv99aNT/GeAl3kJOTvLeNHc1exbOteLu7bgd9e3Z+kVjFelyUCKNBFGsQ5x4tLt/PbBeuJMONP1w3kmsEpmoRCmhUFukg9du4/wk/+uZrFm4v5So/2/P7aAaQkatxyaX4U6CLH4ZxjTlY+v3l9PT7nePCqftw4PF1X5dJsKdBFjqFgXxn3vrqGxZuLGdktiT9cO4D0pHivyxI5IQW6SB0+n+OFJdv4w5sbMODBq/pxw7B0dRKSoKBAF/HLLSzlnn+uJnvbPs7vmcxvr+5HaltdlUvwUKBL2Kus9vHk+1t47N1c4qIj9ASLBC0FuoS15dv2ce+rq9m0p5TLB3Til1f0JTlBz5VLcFKgS1g6cKSKh97cwKxl2+nUOla9PSUkKNAlrDjneG31Lu5/LYe9hyv49rlduWtUT1rF6EdBgp8+xRI28opKue8/6/gwt5j+KW2YPmko/VI0MqKEDgW6hLwjlTX8/b1c/vF+HjFRLbh/TF9uHN6FCD2KKCGmQYFuZqOBvwIRwDPOud8f9f4PgSlANVAEfNs5ty3AtYqcFOccC9ft4YH5OezYf4Srzu7MTy87izMSNIuQhKZ6A93MIoDHgVFAAZBlZvOcczl1NlsJZDrnyszsO8BDwNjGKFikIXILS7l/fg4fbCqiV4cE5kwdwYhuSV6XJdKoGnKFPgzIdc7lAZjZHGAM8EWgO+cW1dl+CTA+kEWKNNTB8ioefXsz0z/eSlxUBL+4vA83j+xCVEQLr0sTaXQNCfQUIL/OcgEw/ATbTwbeONYbZjYVmAqQnp7ewBJF6lfjc7ycnc/D/91IyeFKvjUkjbsv7qVnyiWsBPSmqJmNBzKBrx7rfefcU8BTAJmZmS6Q+5bw9cmWEh6Yn0POroMM6dKWaROHMiA10euyRJpcQwJ9B5BWZznVv+5/mNmFwM+ArzrnKgJTnsjx5RWV8vs3NvDfnD2kJMbxt+sHcfmATuqyL2GrIYGeBfQws67UBvk44Ia6G5jZIOAfwGjnXGHAqxSpY+/hSh59ZzMzl2wjJrIF/3dxLyaf15XYqAivSxPxVL2B7pyrNrPbgYXUPrY4zTm3zszuB7Kdc/OAPwKtgJf9V0fbnXNXNmLdEoaOVNYw7aPPePK9LRyurGbcsHTuurCn2slF/BrUhu6cWwAsOGrdfXVeXxjgukS+UF3j45XlBTzy9ib2HKzgwrM68JPRvejRIcHr0kSaFfUUlWbLOccba3fz8MKN5BUfZlB6Io/dMJihGe28Lk2kWVKgS7PjnOO9TUX86b8bWbvjID3OaMVTNw1hVJ8OuuEpcgIKdGlWPt5SzJ//u4nsbftIbRvHw9cN5OpBKRp3RaQBFOjSLCzJK+GRtzax9LO9dGgdw4NX9eNbmWlER6qHp0hDKdDFM845luTt5dF3NvNJXgnJCTH88oo+XD8sXY8gipwCBbo0OeccH+YW87d3clm2dS/JCTH8/LKzGD+ii4Jc5DQo0KXJ+HyOt9fv4fFFuawqOEDH1rH8+sq+jB2apiAXCQAFujS6qhof8z7dyT8+2MKmPaWkt4vnd9f055rBKcREKshFAkWBLo2mtKKal7LyeXZxHjsPlNOrQwJ/GXs2lw/oRKSGsxUJOAW6BNyeg+U899FWXly6jUPl1QzLaMdvru7PBb2S9Ry5SCNSoEvArCk4wLMf5jF/9S58znFJv05M+UpXBqW39bo0kbCgQJfTUlXjY+G63Tz/8Vaytu6jVUwkN4/MYOI5GaQnxXtdnkhYUaDLKSk6VMGcZdt5cel2dh8sJ61dHD+/7CzGDk0jITbK6/JEwpICXRrMOUfW1n28sGQbb67dRVWN4ys92vPgVf34Wu8z1D1fxGMKdKnX/rJK/rliB7OXbSe3sJSE2EhuGpHBjSPSOTO5ldfliYifAl2OyedzfJJXwktZ+by5bjeV1T7OTkvkoWsHcPnATsRH66Mj0tzop1L+R/7eMl5dsYNXVuSTv/cIrWMjGTc0jXFD0+nTubXX5YnICSjQhUPlVbyxZjevrixgSd5ezGBktyTuvqgXF/ftqG75IkFCgR6mKqpr+GBTMf/5dAdv5eyhotpHRlI8PxrVk6sHp5DaVo8cigQbBXoYqa7xsSRvL6+t2skba3dxsLyatvFRjB2axtWDUjg7LVE9OUWCmAI9xFXV+FiSV8KCNbtYuG4Pew9X0jI6gov7duSKsztzXvf2RGlcFZGQoEAPQWWV1SzeXMzCdbt5O2cPB8uriY+O4MKzOnBp/05c0CtZ7eIiIUiBHiL2HCxn0YZC3l6/h8Wbi6mo9tE6NpIL+3RgdN+OnN9TIS4S6hToQarG5/g0fz/vbyzk3Y2FrN1xEICUxDiuH5bOqD4dGNa1nZpTRMKIAj2I5O8t46PcYhbnFvPh5mIOHKmihcGg9Lb838W9+MZZZ9CrQ4JubIqEKQV6M1Z4qJxPtpSwJK+ET7aUsLWkDIAOrWMY1acDF/RK5rzu7UmMj/a4UhFpDhTozYRzjvy9R8jetpdln9V+5RUfBiAhJpLh3dpx88gMvtKjPd3PaKWrcBH5EgW6R45U1rBmxwE+zd/Hyu37yd62j6JDFQAkxEYyLKMdY4emMfLMJPp2bqORDEWkXgr0JnCksob1uw+ybudB1hTsZ3XBATYXllLjcwCktYvj3DOTGJLRjiHpbenVMUEBLiInTYEeQD6fY8f+I2zcfYiNew6xYfch1u86SF5RKf7spm18FP1TE7nwrA6cnZbI2emJtG8V423hIhISFOin4MCRKraXlPFZyWHyikrJKzrMlqJSthSVUl7l+2K71LZx9O6YwKX9O9G3c2v6dm5NSmKc2r9FpFE0KNDNbDTwVyACeMY59/uj3o8BZgBDgBJgrHNua2BLbRo1PkfJ4Qr2HKhg54Ej7D5Qzs79RyjYd4T8fWVs31vG/rKqL7Y3g85t4uiW3JIR3brQ44xW9OjQip4dEjQVm4g0qXoD3cwigMeBUUABkGVm85xzOXU2mwzsc851N7NxwB+AsY1RcENV1fgoq6ihrKqawxU1HCqv4lB5NYfKq9l/pJL9ZVXsO1zJ3rJKSkorKTlcQeHBCkoOV37Rtv256MgWpCbGkdI2jkv7d6JLu3i6JLUko308GUkt1QNTRJqFhlyhDwNynXN5AGY2BxgD1A30McCv/K9fAR4zM3PO/W8yBsDcrHyeWpyHzzmcq72irq7xUe1zVNX4qKiu/To6lI8lNqoFSS1jSGoVTftWMfTp1JozEmI5o3UMHVvH0qlNHJ0SY2kXH00L3aQUkWauIYGeAuTXWS4Ahh9vG+dctZkdAJKA4robmdlUYCpAenr6KRWcGB/l7w0JLcxoYRAZ0YLIFkZkhBEbGUFMVAtiIyOIi46gZUwk8dERJMRG0jo2ioTYKBLjo2gTF6UraxEJKU16U9Q59xTwFEBmZuYpXb1f1LcjF/XtGNC6RERCQUNGbtoBpNVZTvWvO+Y2ZhYJtKH25qiIiDSRhgR6FtDDzLqaWTQwDph31DbzgAn+198E3m2M9nMRETm+eptc/G3itwMLqX1scZpzbp2Z3Q9kO+fmAc8CL5hZLrCX2tAXEZEm1KA2dOfcAmDBUevuq/O6HLgusKWJiMjJ0OwHIiIhQoEuIhIiFOgiIiFCgS4iEiLMq6cLzawI2HaKf709R/VCDRPheNzheMwQnscdjscMJ3/cXZxzycd6w7NAPx1mlu2cy/S6jqYWjscdjscM4Xnc4XjMENjjVpOLiEiIUKCLiISIYA30p7wuwCPheNzheMwQnscdjscMATzuoGxDFxGRLwvWK3QRETmKAl1EJEQEXaCb2Wgz22hmuWZ2j9f1NAYzSzOzRWaWY2brzOxO//p2ZvaWmW32/9nW61oDzcwizGylmc33L3c1s6X+8/2SfwjnkGJmiWb2ipltMLP1ZjYyTM71Xf7P91ozm21msaF2vs1smpkVmtnaOuuOeW6t1qP+Y19tZoNPdn9BFeh1Jqy+BOgDXG9mfbytqlFUAz9yzvUBRgDf8x/nPcA7zrkewDv+5VBzJ7C+zvIfgEecc92BfdROSB5q/gq86ZzrDQyk9vhD+lybWQpwB5DpnOtH7dDcn08wH0rnezow+qh1xzu3lwA9/F9TgSdOdmdBFejUmbDaOVcJfD5hdUhxzu1yzq3wvz5E7Q94CrXH+rx/s+eBq7ypsHGYWSpwGfCMf9mAr1M78TiE5jG3Ac6ndk4BnHOVzrn9hPi59osE4vyznMUDuwix8+2c+4DaOSLqOt65HQPMcLWWAIlm1ulk9hdsgX6sCatTPKqlSZhZBjAIWAp0cM7t8r+1G+jgUVmN5S/AjwGffzkJ2O+cq/Yvh+L57goUAc/5m5qeMbOWhPi5ds7tAB4GtlMb5AeA5YT++Ybjn9vTzrdgC/SwYmatgH8CP3DOHaz7nn+Kv5B55tTMLgcKnXPLva6liUUCg4EnnHODgMMc1bwSaucawN9uPIba/9A6Ay35ctNEyAv0uQ22QG/IhNUhwcyiqA3zF51zr/pX7/n8V3oKHlsAAAFYSURBVDD/n4Ve1dcIzgWuNLOt1DalfZ3atuVE/6/kEJrnuwAocM4t9S+/Qm3Ah/K5BrgQ+Mw5V+ScqwJepfYzEOrnG45/bk8734It0BsyYXXQ87cdPwusd879uc5bdSfjngD8p6lrayzOuXudc6nOuQxqz+u7zrkbgUXUTjwOIXbMAM653UC+mfXyr/oGkEMIn2u/7cAIM4v3f94/P+6QPt9+xzu384Cb/U+7jAAO1GmaaRjnXFB9AZcCm4AtwM+8rqeRjvE8an8NWw186v+6lNo25XeAzcDbQDuva22k478AmO9/3Q1YBuQCLwMxXtfXCMd7NpDtP9//BtqGw7kGfg1sANYCLwAxoXa+gdnU3iOoova3scnHO7eAUfsU3xZgDbVPAJ3U/tT1X0QkRARbk4uIiByHAl1EJEQo0EVEQoQCXUQkRCjQRURChAJdRCREKNBFRELE/wNXTID1UOAQawAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 283
        },
        "id": "0VJ3yT1oIYul",
        "outputId": "8476f5f2-edd9-4849-944f-adb0c5ebccfc"
      },
      "source": [
        "import matplotlib.pyplot as plt\n",
        "plt.plot(dirichlet_pdf(1, np.linspace(0,1,100)))"
      ],
      "execution_count": 40,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[<matplotlib.lines.Line2D at 0x7f00e7671a90>]"
            ]
          },
          "metadata": {},
          "execution_count": 40
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deXhV5bn+8e+TGULClDCFeZBBFMEg1KnORe2BarWitYpVOO1Pq621p/rTY9W2p9b22MmhB61aFUWttqBSqQO2thVkHsIY5jAGCCFkHp7zx956UiSygZ2s7L3vz3XtK3sN2etZrOTmzbvetZa5OyIiEvuSgi5ARESiQ4EuIhInFOgiInFCgS4iEicU6CIicSIlqA3n5OR43759g9q8iEhMWrhw4R53zz3cssACvW/fvixYsCCozYuIxCQz29zUMnW5iIjECQW6iEicUKCLiMQJBbqISJxQoIuIxIkjBrqZPWVmu81sRRPLzcx+bWaFZrbMzEZFv0wRETmSSFrozwDjPmP5xcCg8GsK8PjxlyUiIkfriIHu7n8D9n3GKhOAZz1kLtDBzLpHq0ARkXjh7vz4zZWs3VXWLJ8fjT70PGBro+mi8LxPMbMpZrbAzBYUFxdHYdMiIrHjrRU7eeKDjSwvKm2Wz2/Rk6LuPtXd8909Pzf3sFeuiojEpfoG57/fXsuA3Ey+NPKwbd7jFo1A3wb0ajTdMzxPRETCZi7dRuHug9x+4WCSk6xZthGNQJ8JXBce7TIWKHX3HVH4XBGRuFBb38Av3l7HsO7ZXDy8W7Nt54g35zKzF4FzgBwzKwJ+AKQCuPtvgVnAJUAhUAHc0FzFiojEolcWFLFlXwVPTconqZla5xBBoLv71UdY7sDNUatIRCSOVNXW85v31jGqdwfOHdylWbelK0VFRJrR83M3s6O0iju+MBiz5mudgwJdRKTZlFXV8uicQs4alMPpA3KafXsKdBGRZvLkBxspqajle18Y3CLbU6CLiDSDvQerefKDDVxyUjdO7tmhRbapQBcRaQaPzllPVV0D372oZVrnoEAXEYm6bfsreX7uZq4Y1ZMBue1abLsKdBGRKHv4L2vB4LYLBrXodhXoIiJRtGrHAV5bXMQNp/elR4c2LbptBbqISBQ99NZqstJT+OY5A1p82wp0EZEo+XD9XuasKeb/nTuQDm3TWnz7CnQRkShwdx58azXd22cw6fS+gdSgQBcRiYJZy3eydOt+vnPhCWSkJgdSgwJdROQ41dQ18NDs1QzumsWXR/UMrA4FuojIcXp+7mY2763gzkuGNNvDKyKhQBcROQ6llbX8+r11nDGwM+ecEOyjNRXoIiLH4bE5hZRW1vL/Lxna7LfHPRIFuojIMSoqqeDpf27ispF5nNijfdDlKNBFRI7VQ2+twYA7WvAGXJ9FgS4icgwWbSlh5tLtTD6rf4tf4t8UBbqIyFFyd370xkpys9L5RgCX+DdFgS4icpTeXL6DRVv2c8dFJ9AuPSXocj6hQBcROQpVtfU8+OfVDOmWxRWn9gq6nH+hQBcROQpP/2MTRSWV3HPpsEAvIjocBbqISIR2l1XxyHvruGBoF84clBN0OZ+iQBcRidDP3lpDTX0Dd186LOhSDkuBLiISgeVFpfxhURE3nNGPfjmZQZdzWAp0EZEjcHfuf72ATm3TuOW8gUGX0yQFuojIEby+bAcLNpfwvS8MJjsjNehymqRAFxH5DBU1dfxk1ipO7JHNlfmta5jioVrPiHgRkVbosTnr2VFaxW+uHtnqhikeSi10EZEmbN5bztS/beBLp/Qgv2+noMs5oogC3czGmdkaMys0szsPs7y3mc0xs8VmtszMLol+qSIiLeuHb6wiJdm465KhQZcSkSMGupklA48CFwPDgKvN7NBBmPcAL7v7SGAi8Fi0CxURaUnvr9nNO6t28a3zBtE1OyPociISSQv9NKDQ3Te4ew0wHZhwyDoOZIfftwe2R69EEZGWVV1Xz/2vr6RfTiZfP7Nv0OVELJJAzwO2NpouCs9r7D7gWjMrAmYB3zrcB5nZFDNbYGYLiouLj6FcEZHm9+QHG9m4p5wf/Nsw0lOSgy4nYtE6KXo18Iy79wQuAZ4zs099trtPdfd8d8/PzQ32YaoiIoezbX8lv3lvHeNO7MY5g7sEXc5RiSTQtwGNB1/2DM9r7EbgZQB3/xDIAFrfnWtERI7gh6+vBOA//6113q/ls0QS6POBQWbWz8zSCJ30nHnIOluA8wHMbCihQFefiojElL+uLeatgp1867xB5LWSx8odjSMGurvXAbcAs4FVhEazFJjZA2Y2Przad4HJZrYUeBGY5O7eXEWLiERbVW09P5ixgn45mdx0Vr+gyzkmEV0p6u6zCJ3sbDzv3kbvVwJnRLc0EZGW89u/rmfT3gqev3FMTJ0IbUxXiopIwtu4p5zH5qxn/IgerfLBFZFSoItIQnN37p2xgvSUJO75YmxcEdoUBbqIJLQ3lu3gg3V7uOMLg+mSFRtXhDZFgS4iCau0spYH3ljJSXntuXZsn6DLOW66fa6IJKyH3lrN3oPVPD1pdKu/NW4k1EIXkYS0cPM+ps3bwg1n9GN4Xvugy4kKBbqIJJyaugbuem05PdpncPuFJwRdTtSoy0VEEs4TH2xg7a6DPHldPpnp8RODaqGLSELZuKecX727jouHd+OCYV2DLieqFOgikjAaGpw7X11GekoS948/Mehyok6BLiIJ46UFW5m3cR93XzKULjHyFKKjoUAXkYSw60AV/zVrFZ/r35mrRvc68jfEIAW6iMS9jy/vr6lr4CeXn4RZ7I85PxwFuojEvVnLdzK7YBffufAE+uZkBl1Os1Ggi0hc23uwmntnrODknu256czYvM95pOJnAKaIyGHc//pKDlTVMu2KMaQkx3cbNr73TkQS2l8KdjJz6XZuOXcQQ7plB11Os1Ogi0hc2l9Rwz1/WsGQbll885wBQZfTItTlIiJx6f7XV7KvvIanJo0mLSUx2q6JsZciklBmF+zkj4u3cfO5A+PmToqRUKCLSFzZV17D3X9czrDu2dx87sCgy2lR6nIRkbjyg5kFlFbW8tyNYxKmq+VjibW3IhLXXl+6ndeXbufW8wYxtHv8j2o5lAJdROLCrgNV/OeMFYzo1SFhRrUcSoEuIjHP3fn+q8uoqq3n4a+MiPsLiJqSmHstInHlxY+28v6aYu66eCgDctsFXU5gFOgiEtM27SnnR2+u5MyBOXxtbJ+gywmUAl1EYlZdfQPffmkJKUnGQ1ecTFJSfN4WN1IatigiMes37xWyZOt+HrlmJD06tAm6nMCphS4iMWnh5hIemVPI5SPz+OLJPYIup1WIKNDNbJyZrTGzQjO7s4l1vmJmK82swMxeiG6ZIiL/52B1Hbe/vIRu2RncNyH+HvZ8rI7Y5WJmycCjwIVAETDfzGa6+8pG6wwC7gLOcPcSM+vSXAWLiPxgRgFb91UwfcrnyM5IDbqcViOSFvppQKG7b3D3GmA6MOGQdSYDj7p7CYC7745umSIiITOWbOPVRUXccu5ATuvXKehyWpVIAj0P2Npouig8r7ETgBPM7B9mNtfMxh3ug8xsipktMLMFxcXFx1axiCSsrfsquOePKzi1T0duPX9Q0OW0OtE6KZoCDALOAa4GnjCzDoeu5O5T3T3f3fNzc3OjtGkRSQS19Q3cOn0xGPzyqlMS9mrQzxLJv8g2oFej6Z7heY0VATPdvdbdNwJrCQW8iEhU/OLttSzesp8fX3YSvTq1DbqcVimSQJ8PDDKzfmaWBkwEZh6yzp8Itc4xsxxCXTAboliniCSwD9YV8/hf13NVfi/Gj9AQxaYcMdDdvQ64BZgNrAJedvcCM3vAzMaHV5sN7DWzlcAc4Hvuvre5ihaRxLG7rIrvvLSEgbntuG+8hih+loiuFHX3WcCsQ+bd2+i9A7eHXyIiUdHQ4Nz+0lIOVtfxwuSxtElLDrqkVk2X/otIq/XY+4X8vXAPD15+Eid0zQq6nFZPp4lFpFX65/o9PPz2Wiac0oOrRvc68jeIAl1EWp/dZVXc+uIS+uVk8l+XnYRZYt9FMVLqchGRVqW+wbn1xcUcrK5l2k1jyExXTEVK/1Ii0qo8/PYa5m7Yx8+vHMHgbuo3PxrqchGRVuPtlbt4dE5ovPkVp/YMupyYo0AXkVZh055ybn95CcPzsrlft8Q9Jgp0EQlcZU0933h+IUlmPP7VU8lI1XjzY6E+dBEJlLtz9x+Xs2ZXGU9NGq37tBwHtdBFJFDP/HMTry3exm3nD+LcwXo2zvFQoItIYOZt2MuP3lzFBUO7cut5ukHr8VKgi0ggdpRWcvMLi+jTqS0PXzWCpCRdPHS81IcuIi2uqraebzy3kMqaeqZPGavngkaJAl1EWpS7c+ery1haVMrUr53KwC66eCha1OUiIi1q6t828Kcl2/nuhSdw0Yndgi4nrijQRaTFzFmzmwffWs2lJ3XnlvMGBl1O3FGgi0iLWLurjG+9sJih3bL52ZUn6w6KzUCBLiLNbl95DTf+fj4Zqck8eX0+bdN0+q456F9VRJpVTV0D33huIbsOVPPSlLH06NAm6JLillroItJsPr6s/6NNodvhjuzdMeiS4poCXUSazWPvr+eVhUXcev4gxo/oEXQ5cU+BLiLNYsaSbfxs9houG5nHdy7QZf0tQYEuIlE3f9M+vvfKMk7r24kHv6xngrYUBbqIRNX64oNMfnYBeR3b8D9fO5X0FN3bvKUo0EUkaorLqpn09Eckm/HMDaPpmJkWdEkJRcMWRSQqyqvr+Poz8ykuq2b6lM/Rp3Nm0CUlHLXQReS41dY3cMsLiyjYXsojV4/ilF4dgi4pIamFLiLHJXT3xOXMWVPMjy8bzgXDugZdUsJSC11EjstDs9fw6qIivn3BIL46pk/Q5SQ0BbqIHLPf/X0jj7+/nq+O6c1t52usedAU6CJyTF5dWMQP31jJuBO78cCE4Rpr3gpEFOhmNs7M1phZoZnd+RnrfdnM3Mzyo1eiiLQ2fynYyX+8uozTB3TmlxNPIVnPA20VjhjoZpYMPApcDAwDrjazYYdZLwu4DZgX7SJFpPX45/o93PLiYobntWfqdflkpOrCodYikhb6aUChu29w9xpgOjDhMOv9EPgpUBXF+kSkFVm0pYTJv19A385teWbSaNqla6BcaxJJoOcBWxtNF4XnfcLMRgG93P3Nz/ogM5tiZgvMbEFxcfFRFysiwSnYXsqkpz4iJyud524co6tAW6HjPilqZknAw8B3j7Suu09193x3z8/NzT3eTYtICyncXcbXfvcR7dJTmHbTGLpmZwRdkhxGJIG+DejVaLpneN7HsoDhwPtmtgkYC8zUiVGR+LBxTznXPDGP5CRj2uSx9OzYNuiSpAmRBPp8YJCZ9TOzNGAiMPPjhe5e6u457t7X3fsCc4Hx7r6gWSoWkRazeW85V0+dS12DM+2mMfTL0f1ZWrMjBrq71wG3ALOBVcDL7l5gZg+Y2fjmLlBEgrF1XwXXPDGPqrp6pt00hhO6ZgVdkhxBRKeo3X0WMOuQefc2se45x1+WiASpqKSCq5+Yy8HqOqbdNIah3bODLkkioDFHIvIvtu6rYOLUuZRV1fL8TWMYntc+6JIkQgp0EfnEx33m5TX1vDB5rMI8xuheLiIChEazTJw6l4raUJ+5wjz2qIUuIqzdVcZXn5xHfYPzwk1jGdZDfeaxSC10kQS3YlspV/3Phxjw0hSFeSxTC10kgS3cXMKkpz8iOyOVaTeNoa/Gmcc0tdBFEtQH64q59sl5dM5M46V/H6swjwNqoYskoD8v38Gt0xczILcdz954Gl2ydG+WeKBAF0kwL83fwl2vLWdk7448NWk07dukBl2SRIkCXSRBuDuPzink539Zy9kn5PLba0fRNk0REE90NEUSQEODc//rBfz+w81cNjKPh644mdRknUKLNwp0kThXVVvP7S8vYdbynUw+qx93XTyUJD0DNC4p0EXiWEl5DZOfXcCCzSXcc+lQbjqrf9AlSTNSoIvEqa37Krj+6Y8oKqnk0WtGcenJ3YMuSZqZAl0kDi3cvI8pzy785MEUo/t2CrokaQEKdJE4M3Ppdu54ZSk92mfw1KTR9M9tF3RJ0kIU6CJxwt359buF/OKdtZzWtxO//dqpdMpMC7osaUEKdJE4UFlTzx1/WMqby3Zw+cg8fvLlk0hPSQ66LGlhCnSRGLejtJLJzy6gYPsB7rp4CFPO7o+ZhiUmIgW6SAybv2kf33x+EZU1dTx5XT7nD+0adEkSIAW6SAxyd6bN28J9Mwvo2bENL0wewwlds4IuSwKmQBeJMVW19dw3s4Dp87dyzuBcfjVxpG6wJYACXSSmbN1XwTenLWTFtgPcfO4Abr9wMMm6jF/CFOgiMWLOmt18e/oSGtx54rp8Lhym/nL5Vwp0kVaurr6Bh99ey2Pvr2dItyx+e+2perqQHJYCXaQV21laxbdeXMT8TSVMHN2L+8afSEaqxpfL4SnQRVqp91bv4o5XllFVW88vrhrBZSN7Bl2StHIKdJFWprqungf/vJqn/7GJod2z+c3VIxnYRfdjkSNToIu0Iut2lXHb9CWs3HGASaf35c6Lh6iLRSKmQBdpBdydZz/czH/NWkVmeopGscgxiSjQzWwc8CsgGXjS3R88ZPntwE1AHVAMfN3dN0e5VpG4tOtAFd9/dRnvrynmnMG5PHTFyXTJygi6LIlBRwx0M0sGHgUuBIqA+WY2091XNlptMZDv7hVm9k3gIeCq5ihYJF64OzOXbufeGQVU19XzwIQT+drYPrqxlhyzSFropwGF7r4BwMymAxOATwLd3ec0Wn8ucG00ixSJN3sOVnPvjBXMWr6TUb078N9fOYV+GlsuxymSQM8DtjaaLgLGfMb6NwJ/PtwCM5sCTAHo3bt3hCWKxI+PW+X3zSygvLqe748L3e5Wl+9LNET1pKiZXQvkA58/3HJ3nwpMBcjPz/dobluktdtZWsU9f1rBO6t2MbJ3B352xckM7KI7JEr0RBLo24BejaZ7huf9CzO7ALgb+Ly7V0enPJHY19DgTJu3mZ++tYba+gbuvmQoXz+zn1rlEnWRBPp8YJCZ9SMU5BOBaxqvYGYjgf8Bxrn77qhXKRKjVu88wN1/XMHCzSWcOTCHH182nD6d1VcuzeOIge7udWZ2CzCb0LDFp9y9wMweABa4+0zgZ0A74JXwGfot7j6+GesWadXKq+v41bvr+N3fN5KdkcJ/XzmCy0flaQSLNKuI+tDdfRYw65B59zZ6f0GU6xKJSe7OrOU7+dGbK9lRWsXE0b34/rghdMxMC7o0SQC6UlQkStbuKuMHMwr4cMNehnbP5pFrRnFqn45BlyUJRIEucpz2V9Twy3fW8dzczbRLT+GHE07kmjF9dNJTWpwCXeQY1dY3MG3uZn7xzjrKqmq5+rTefPeiwXRS94oERIEucpTcndkFO/npW2vYuKecMwfmcM8XhzKkW3bQpUmCU6CLHIX5m/bxk1mrWLRlPwO7tON31+dz3pAuGr0irYICXSQCBdtL+fnsNcxZU0xuVjo/ufwkrjy1JynJSUGXJvIJBbrIZyjcXcYv31nHG8t2kJ2RwvfHDeH60/vQNk2/OtL66KdS5DDWFx/kN++uY8bS7bRJTebmcwcw5ewBtG+TGnRpIk1SoIs0smZnGY/MKeTNZdtJT0lmytn9mXJWfzq3Sw+6NJEjUqCLAEu27ufx9wuZXbCLzLRkJp/dn8ln9SdHQS4xRIEuCcvd+du6Pfz2/fV8uGEv2Rkp3HreQG44o58u1ZeYpECXhFNdV8+Mxdt58u8bWLvrIN2yM7jn0qFMPK037dL1KyGxSz+9kjB2H6hi2rwtTJu3hT0HqxnSLYufXzmC8SN6kJai4YcS+xToEtfcnUVbSnj2w83MWr6D2nrnvCFduOGMvpw5MEcXBElcUaBLXCqrqmXGku08P3czq3eWkZWewlfH9OH60/vqYcwStxToEjdCrfH9TP9oC28s20FlbT3Dumfzk8tPYvyIHmSqf1zinH7CJebtOlDFa4u28YeFW1lfXE7btGQmnNKDq0b34pReHdStIglDgS4x6WB1HbNX7ORPS7bxj8I9NDiM7tuRyWf154sjemi0iiQk/dRLzKiqref9NcW8vmw7767aRVVtA706teHmcwdy+aie6huXhKdAl1bt4xB/a8UO3lm1m4PVdXTKTOPLo3py2cg8Tu3TUV0qImEKdGl19lfU8N7q3by9chd/XVtMRU09HdqmcslJ3fi3ET34XP/Oum2tyGEo0CVw7s764oO8u2o3763ezYLNJdQ3OF2y0rlsZB4XD+/OmP6dSFWIi3wmBboEoqyqlg/X7+Wva4v569piikoqARjSLYt/P7s/F53YjZPz2pOkBy2LREyBLi2iqraeJVv388/CPfy9cA9Li0qpb3Ay05I5fWAO3/j8AM4d0oW8Dm2CLlUkZinQpVlU1NSxeMt+Ptq4j3kb97Joy35q6hpIMji5Zwe++fkBnDEwh1P7dNR9VESiRIEux83d2ba/kkVb9rNocwmLt5SwYvsB6hscMxjWPZvrxvZhbP/OjO7XSU/9EWkmCnQ5arsOVLFiWynLt5WyrKiUpVv3s7e8BoA2qcmc3LM93/h8f0b37cSoPh3JzlCAi7QEBbo0qaaugY17ylm98wCrdpSxcscBVu04QHFZNQBmMDC3HecO6cKInu0Z2bsjQ7plaUihSEAU6EJFTR2b9lSwvvgghbsPUlh8kMJdB1lffJC6BgcgNdkY1CWLswflMjwvm+F57RnaPVuX2Iu0IvptTADuTmllLVv3VbK1pIIt+yrYvLecTXsq2LS3nB2lVZ+sawa9O7VlYG47zhvahSHdsjihaxYDctvp5KVIK6dAj3HuzoHKOnaVVbHrQBU7S0OvHQeq2LG/km37K9lWUkl5Tf2/fF/nzDT6dG7L2P6d6Z+TSf/cdvTLyaR/biYZqckB7Y2IHI+IAt3MxgG/ApKBJ939wUOWpwPPAqcCe4Gr3H1TdEtNDLX1DZRW1lJaWcv+ilr2V9RQUlFLSXkN+ypq2Huwmr0Ha9hTXsOesmqKD1ZTU9fwqc/pnJlGt/YZ9O2cyekDcujZsQ29OrWlV8e29OrUhiydqBSJO0cMdDNLBh4FLgSKgPlmNtPdVzZa7UagxN0HmtlE4KfAVc1RcBDcnfoGpy78qq1roLahgdp6p6augdr6BmrqGqiuq6eq9v++VtbUU1VXT2VN6FVRW095dR3l1eGvNXWUVdVxsLqOsqpaDlTWUVlb32QdKUlG53ZpdMpMJ6ddGgNyM8ltl05uVjpdszPoEv7arX2GWtkiCSiSFvppQKG7bwAws+nABKBxoE8A7gu//wPwiJmZu3sUawXg5flbeeKDDZ9MN95A4815ozceXhb6Co7TEG7UNriHX9DQ4NSHw7shHN4N7tTWR2c30pKTaJueTGZaCm3TkmmXkUJWRgrd22eQnZFKdpsUsjJS6dA2lfZtUsluk0rHtml0bJtKh7ZpZGek6M6CItKkSAI9D9jaaLoIGNPUOu5eZ2alQGdgT+OVzGwKMAWgd+/ex1Rwh7apDOra7l/mGdZ44lNvzUJrmIXmJYXfGEZyUmjaLPQ1Ock++ZqSZCQlGalJRnJSEinJRmqykZqcFH4ZaSlJpCUnk5aSRPrHr9Rk2qQmk5GaREZqMm3SkmmbmqzhfCLSrFr0pKi7TwWmAuTn5x9Ts/eiE7tx0YndolqXiEg8iKTJuA3o1Wi6Z3jeYdcxsxSgPaGToyIi0kIiCfT5wCAz62dmacBEYOYh68wErg+/vwJ4rzn6z0VEpGlH7HIJ94nfAswmNGzxKXcvMLMHgAXuPhP4HfCcmRUC+wiFvoiItKCI+tDdfRYw65B59zZ6XwVcGd3SRETkaGjYhYhInFCgi4jECQW6iEicUKCLiMQJC2p0oZkVA5uP8dtzOOQq1ASRiPudiPsMibnfibjPcPT73cfdcw+3ILBAPx5mtsDd84Ouo6Ul4n4n4j5DYu53Iu4zRHe/1eUiIhInFOgiInEiVgN9atAFBCQR9zsR9xkSc78TcZ8hivsdk33oIiLyabHaQhcRkUMo0EVE4kTMBbqZjTOzNWZWaGZ3Bl1PczCzXmY2x8xWmlmBmd0Wnt/JzN42s3Xhrx2DrjXazCzZzBab2Rvh6X5mNi98vF8K38I5rphZBzP7g5mtNrNVZva5BDnW3wn/fK8wsxfNLCPejreZPWVmu81sRaN5hz22FvLr8L4vM7NRR7u9mAr0Rg+svhgYBlxtZsOCrapZ1AHfdfdhwFjg5vB+3gm86+6DgHfD0/HmNmBVo+mfAr9w94FACaEHksebXwFvufsQYASh/Y/rY21mecCtQL67Dyd0a+6PHzAfT8f7GWDcIfOaOrYXA4PCrynA40e7sZgKdBo9sNrda4CPH1gdV9x9h7svCr8vI/QLnkdoX38fXu33wJeCqbB5mFlP4FLgyfC0AecRevA4xOc+twfOJvRMAdy9xt33E+fHOiwFaBN+yllbYAdxdrzd/W+EnhHRWFPHdgLwrIfMBTqYWfej2V6sBfrhHlidF1AtLcLM+gIjgXlAV3ffEV60E+gaUFnN5ZfAfwAN4enOwH53rwtPx+Px7gcUA0+Hu5qeNLNM4vxYu/s24OfAFkJBXgosJP6PNzR9bI8732It0BOKmbUDXgW+7e4HGi8LP+IvbsacmtkXgd3uvjDoWlpYCjAKeNzdRwLlHNK9Em/HGiDcbzyB0H9oPYBMPt01EfeifWxjLdAjeWB1XDCzVEJhPs3dXwvP3vXxn2Dhr7uDqq8ZnAGMN7NNhLrSziPUt9wh/Cc5xOfxLgKK3H1eePoPhAI+no81wAXARncvdvda4DVCPwPxfryh6WN73PkWa4EeyQOrY1647/h3wCp3f7jRosYP474emNHStTUXd7/L3Xu6e19Cx/U9d/8qMIfQg8chzvYZwN13AlvNbHB41vnASuL4WIdtAcaaWdvwz/vH+x3XxzusqWM7E7guPNplLFDaqGsmMu4eUy/gEmAtsB64O+h6mmkfzyT0Z9gyYEn4dQmhPuoplBAAAACJSURBVOV3gXXAO0CnoGttpv0/B3gj/L4/8BFQCLwCpAddXzPs7ynAgvDx/hPQMRGONXA/sBpYATwHpMfb8QZeJHSOoJbQX2M3NnVsASM0im89sJzQCKCj2p4u/RcRiROx1uUiIiJNUKCLiMQJBbqISJxQoIuIxAkFuohInFCgi4jECQW6iEic+F/yOUMvRgy6SwAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "plTTdeiz08ac"
      },
      "source": [
        "Multivariate normal distribution for continuous variables x1 and x2"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "KIXs1GXN23JB",
        "outputId": "f406f8fb-8689-4205-eaf0-2f0d25e1417f"
      },
      "source": [
        "from scipy.stats import multivariate_normal\n",
        "\n",
        "std_normal_pdf = lambda *x: multivariate_normal.pdf(x, [0, 0], [[1, 0], [0, 1]])\n",
        "std_normal = ContinuousFactor(['x1', 'x2'], std_normal_pdf)\n",
        "std_normal.scope(), std_normal.assignment([1, 1])"
      ],
      "execution_count": 41,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "(['x1', 'x2'], 0.05854983152431917)"
            ]
          },
          "metadata": {},
          "execution_count": 41
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "k8DKJ_NY1C7n"
      },
      "source": [
        "Most common distribution for continuous variables is Linear Gaussian CPD"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "kJf7MZM01a18",
        "outputId": "0012cb76-105e-4c71-d949-9b9026749426"
      },
      "source": [
        "from pgmpy.factors.continuous import LinearGaussianCPD\n",
        "cpd = LinearGaussianCPD('Y',  [9, -2, 3, 8], 9.6, ['X1', 'X2', 'X3'])\n",
        "print(cpd)"
      ],
      "execution_count": 9,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "P(Y | X1, X2, X3) = N(9*X1 + -2*X2 + 3*X3 + 9; 9.6)\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "fQbzdwLJ1n-3",
        "outputId": "1ad618d7-02e9-49e4-c7da-799fe9c57c9e"
      },
      "source": [
        "cpd.variable, cpd.evidence, cpd.mean"
      ],
      "execution_count": 17,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "('Y', ['X1', 'X2', 'X3'], [9, -2, 3, 8])"
            ]
          },
          "metadata": {},
          "execution_count": 17
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nnhUCM5R3dkl"
      },
      "source": [
        "Joint Linear Gaussian CPDs"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "D0SiJwBE3T6N",
        "outputId": "957625e8-b0f9-4cda-8803-f8118dd94011"
      },
      "source": [
        "from pgmpy.models import LinearGaussianBayesianNetwork\n",
        "\n",
        "model = LinearGaussianBayesianNetwork([('x1', 'x2'), ('x2', 'x3')])\n",
        "cpd1 = LinearGaussianCPD('x1', [1], 4)\n",
        "cpd2 = LinearGaussianCPD('x2', [-5, 0.5], 4, ['x1'])\n",
        "cpd3 = LinearGaussianCPD('x3', [4, -1], 3, ['x2'])\n",
        "# This is a hack due to a bug in pgmpy (LinearGaussianCPD\n",
        "# doesn't have `variables` attribute but `add_cpds` function\n",
        "# wants to check that...)\n",
        "cpd1.variables = [*cpd1.evidence, cpd1.variable]\n",
        "cpd2.variables = [*cpd2.evidence, cpd2.variable]\n",
        "cpd3.variables = [*cpd3.evidence, cpd3.variable]\n",
        "model.add_cpds(cpd1, cpd2, cpd3)\n",
        "jgd = model.to_joint_gaussian()\n",
        "jgd.variables"
      ],
      "execution_count": 18,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "['x1', 'x2', 'x3']"
            ]
          },
          "metadata": {},
          "execution_count": 18
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "ct0gB1N03lTY",
        "outputId": "81addcc4-884a-4a5f-8fb2-27d5933df676"
      },
      "source": [
        "jgd.mean"
      ],
      "execution_count": 19,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "array([[  1.],\n",
              "       [-10.],\n",
              "       [-36.]])"
            ]
          },
          "metadata": {},
          "execution_count": 19
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "g2ctNccQ3weB",
        "outputId": "30f3f616-7de6-405f-a937-9c77e298ad65"
      },
      "source": [
        "cpd2.pdf((-5,0.2))"
      ],
      "execution_count": 30,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "array([0.19947114, 0.00679148])"
            ]
          },
          "metadata": {},
          "execution_count": 30
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "MzXOhhOHxfEZ",
        "outputId": "7eeb0322-95f2-4bf7-d7f5-711f34997584"
      },
      "source": [
        "from pgmpy.models import BayesianModel\n",
        "from pgmpy.factors.discrete import TabularCPD\n",
        "\n",
        "# Defining the model structure. We can define the network by just passing a list of edges.\n",
        "model = BayesianModel([('D', 'G'), ('I', 'G'), ('G', 'L'), ('I', 'S')])"
      ],
      "execution_count": 2,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "/usr/local/lib/python3.7/dist-packages/pgmpy/models/BayesianModel.py:10: FutureWarning: BayesianModel has been renamed to BayesianNetwork. Please use BayesianNewtork class, BayesianModel will be removed in future.\n",
            "  FutureWarning,\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "P1SJnHl-zLvn",
        "outputId": "d1e0dd54-1d1d-4b2b-d4d1-d60e820327a3"
      },
      "source": [
        "model.check_model()"
      ],
      "execution_count": 11,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "True"
            ]
          },
          "metadata": {},
          "execution_count": 11
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "vnPKGzbPxuoC"
      },
      "source": [
        "\n",
        "# Defining individual CPDs.\n",
        "cpd_d = TabularCPD(variable='D', variable_card=2, values=[[0.6], [0.4]])\n",
        "cpd_i = TabularCPD(variable='I', variable_card=2, values=[[0.7], [0.3]])\n",
        "\n",
        "# The representation of CPD in pgmpy is a bit different than the CPD shown in the above picture. In pgmpy the colums\n",
        "# are the evidences and rows are the states of the variable. So the grade CPD is represented like this:\n",
        "#\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | diff    | intel_0 | intel_0 | intel_1 | intel_1 |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | intel   | diff_0  | diff_1  | diff_0  | diff_1  |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_0 | 0.3     | 0.05    | 0.9     | 0.5     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_1 | 0.4     | 0.25    | 0.08    | 0.3     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "#    | grade_2 | 0.3     | 0.7     | 0.02    | 0.2     |\n",
        "#    +---------+---------+---------+---------+---------+\n",
        "\n",
        "cpd_g = TabularCPD(variable='G', variable_card=3,\n",
        "                   values=[[0.3, 0.05, 0.9,  0.5],\n",
        "                           [0.4, 0.25, 0.08, 0.3],\n",
        "                           [0.3, 0.7,  0.02, 0.2]],\n",
        "                  evidence=['I', 'D'],\n",
        "                  evidence_card=[2, 2])\n",
        "\n",
        "cpd_l = TabularCPD(variable='L', variable_card=2,\n",
        "                   values=[[0.1, 0.4, 0.99],\n",
        "                           [0.9, 0.6, 0.01]],\n",
        "                   evidence=['G'],\n",
        "                   evidence_card=[3])\n",
        "\n",
        "cpd_s = TabularCPD(variable='S', variable_card=2,\n",
        "                   values=[[0.95, 0.2],\n",
        "                           [0.05, 0.8]],\n",
        "                   evidence=['I'],\n",
        "                   evidence_card=[2])\n",
        "\n",
        "# Associating the CPDs with the network\n",
        "model.add_cpds(cpd_d, cpd_i, cpd_g, cpd_l, cpd_s)\n",
        "\n",
        "# check_model checks for the network structure and CPDs and verifies that the CPDs are correctly\n",
        "# defined and sum to 1.\n",
        "model.check_model()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "UDKj_fcS1YFS"
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}