diff --git a/.github/workflows/publish-to-pypi.yml b/.github/workflows/publish-to-pypi.yml new file mode 100644 index 0000000..f74e567 --- /dev/null +++ b/.github/workflows/publish-to-pypi.yml @@ -0,0 +1,51 @@ +name: Publish Python ๐Ÿ distribution ๐Ÿ“ฆ to PyPI and TestPyPI + +on: push + +jobs: + build: + name: Build distribution ๐Ÿ“ฆ + runs-on: ubuntu-latest + + steps: + - uses: actions/checkout@v4 + - name: Set up Python + uses: actions/setup-python@v4 + with: + python-version: "3.x" + + - name: Install pypa/build + run: >- + python3 -m + pip install + build + --user + - name: Build a binary wheel and a source tarball + run: python3 -m build + - name: Store the distribution packages + uses: actions/upload-artifact@v3 + with: + name: python-package-distributions + path: dist/ + + publish-to-pypi: + name: >- + Publish Python ๐Ÿ distribution ๐Ÿ“ฆ to PyPI + if: startsWith(github.ref, 'refs/tags/') # only publish to PyPI on tag pushes + needs: + - build + runs-on: ubuntu-latest + environment: + name: pypi + url: https://pypi.org/p/random_events # Replace with your PyPI project name + permissions: + id-token: write # IMPORTANT: mandatory for trusted publishing + + steps: + - name: Download all the dists + uses: actions/download-artifact@v3 + with: + name: python-package-distributions + path: dist/ + - name: Publish distribution ๐Ÿ“ฆ to PyPI + uses: pypa/gh-action-pypi-publish@release/v1 diff --git a/.github/workflows/test-and-build.yml b/.github/workflows/test-and-build.yml new file mode 100644 index 0000000..f192353 --- /dev/null +++ b/.github/workflows/test-and-build.yml @@ -0,0 +1,81 @@ +name: Run tests and build the Python distribution + +on: + push: + branches: + - master + + # -------------------------------------------------------------------------------------------------------------------- + + workflow_call: + inputs: + version: + required: true + type: string + + python-versions: + required: true + type: string + default: "cp38-cp38" + +# ---------------------------------------------------------------------------------------------------------------------- + +defaults: + run: + shell: bash + working-directory: . + +# ---------------------------------------------------------------------------------------------------------------------- + +concurrency: + group: 'random-events-building-and-deployment' + cancel-in-progress: true + +# ---------------------------------------------------------------------------------------------------------------------- + +jobs: + test-and-build: + name: Build and Test the Python distribution + runs-on: ubuntu-22.04 + + steps: + + - name: Checkout ๐Ÿ›Ž + uses: actions/checkout@v3 + + # ---------------------------------------------------------------------------------------------------------------- + + - name: Setup python ๐Ÿ + uses: actions/setup-python@v4 + with: + python-version: 3.8 + cache: pip + + # ---------------------------------------------------------------------------------------------------------------- + + - name: Install user dependencies ๐Ÿผ + uses: py-actions/py-dependency-install@v4 + with: + path: "requirements.txt" + + # ---------------------------------------------------------------------------------------------------------------- + + - name: Run Tests ๐ŸŽ“ + run: | + cd test + PYTHONPATH=../src python -m unittest discover + + # ---------------------------------------------------------------------------------------------------------------- + + - name: Install Sphinx dependencies ๐Ÿ“š + uses: py-actions/py-dependency-install@v4 + with: + path: "doc/requirements.txt" + + # ---------------------------------------------------------------------------------------------------------------- + + - name: Build Sphinx documentation ๐Ÿ“ + working-directory: ./doc + run: | + sudo apt install pandoc + make html \ No newline at end of file diff --git a/doc/notebooks/quickstart.ipynb b/doc/notebooks/quickstart.ipynb deleted file mode 120000 index 89ef87d..0000000 --- a/doc/notebooks/quickstart.ipynb +++ /dev/null @@ -1 +0,0 @@ -../../example/quickstart.ipynb \ No newline at end of file diff --git a/doc/notebooks/quickstart.ipynb b/doc/notebooks/quickstart.ipynb new file mode 100644 index 0000000..1a85469 --- /dev/null +++ b/doc/notebooks/quickstart.ipynb @@ -0,0 +1,241 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "source": [ + "# Quickstart to fglib2\n", + "\n", + "First, let us declare four variables with different domains each." + ], + "metadata": { + "collapsed": false + }, + "id": "155b96b89df61810" + }, + { + "cell_type": "code", + "execution_count": 7, + "outputs": [], + "source": [ + "from random_events.variables import Symbolic\n", + "\n", + "x1 = Symbolic('x1', domain=range(2))\n", + "x2 = Symbolic('x2', domain=range(3))\n", + "x3 = Symbolic('x3', domain=range(4))\n", + "x4 = Symbolic('x4', domain=range(5))" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2023-10-31T09:16:33.220472339Z", + "start_time": "2023-10-31T09:16:33.211022186Z" + } + }, + "id": "b857b83f5ae8482c" + }, + { + "cell_type": "markdown", + "source": [ + "Next, let's create random factors for some pairs of variables. We can shortcut the creation of factors by using the `*` operator instead of adding the nodes and edges manually. " + ], + "metadata": { + "collapsed": false + }, + "id": "7ffcf51c1e699fd5" + }, + { + "cell_type": "code", + "execution_count": 8, + "outputs": [], + "source": [ + "from fglib2.graphs import FactorNode\n", + "from fglib2.distributions import Multinomial\n", + "import numpy as np\n", + "\n", + "np.random.seed(420)\n", + "\n", + "f_x1_x2 = FactorNode(Multinomial([x1, x2], np.random.rand(len(x1.domain), len(x2.domain))))\n", + "f_x2_x3 = FactorNode(Multinomial([x2, x3], np.random.rand(len(x2.domain), len(x3.domain))))\n", + "f_x2_x4 = FactorNode(Multinomial([x2, x4], np.random.rand(len(x2.domain), len(x4.domain))))\n", + "\n", + "graph = f_x1_x2 * f_x2_x3 * f_x2_x4" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2023-10-31T09:16:33.258978788Z", + "start_time": "2023-10-31T09:16:33.258596129Z" + } + }, + "id": "97e60dbee16c1dbd" + }, + { + "cell_type": "markdown", + "source": [ + "We can now draw the graph using ordinary networkx functions." + ], + "metadata": { + "collapsed": false + }, + "id": "3efffd0e9cf4903f" + }, + { + "cell_type": "code", + "execution_count": 9, + "outputs": [ + { + "data": { + "text/plain": "
", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAApQAAAHzCAYAAACe1o1DAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/OQEPoAAAACXBIWXMAAA9hAAAPYQGoP6dpAABRkUlEQVR4nO3deVxU9f7H8ffMsAmJAoptkjsqWVmaS1q5ttg1d+2qmZpLVq5Y3mz1/jLrWll6LfW6lLlvaYu5l6apaO47KUKiIouIgMAsvz+8cUXcYAYODK/n4+HjkXPOnHljKm8/55zvMTkcDocAAACAfDIbHQAAAADFG4USAAAATqFQAgAAwCkUSgAAADiFQgkAAACnUCgBAADgFAolAAAAnEKhBAAAgFMolAAAAHAKhRIAAABOoVACAADAKRRKAAAAOIVCCQAAAKdQKAEAAOAUCiUAAACcQqEEAACAUyiUAAAAcAqFEgAAAE6hUAIAAMApFEoAAAA4hUIJAAAAp1AoAQAA4BQKJQAAAJxCoQQAAIBTKJQAAABwCoUSAAAATqFQAgAAwCkUSgAAADiFQgkAAACnUCgBAADgFAolAAAAnEKhBAAAgFMolAAAAHAKhRIAAABOoVACAADAKRRKAAAAOIVCCQAAAKd4GB0AAACguEnNsCoqIVWZVru8PMyqFOQnP++SW6tK7lcOAACQB8fOpmjOtmhtOBKn6MQ0Oa7YZpIUEuirZqHB6t4gRNUrlDYqpiFMDofDcfPdAAAASqaYxDS9sWyfNkXGy2I2yWa/fnX6a3vTauU0tn0dVQz0LcSkxqFQAgAAXMf8iGi9s+KArHbHDYvk1SxmkzzMJr3XNkzd6ocUYMKigUIJAABwDZM2HNP41UedPk546xp6pVl1FyQqurjLGwAA4CrzI6JdUiYlafzqo1oQEe2SYxVVTCgBAACuEJOYppaf/qIMq/2W9rdeTFTKjhXKiD2izDORcmSmq8JzY+Vzz33Z+3h7mLV22GNue00lE0oAAIArvLFsn6x5uF7SmvCnLmxdLFtKgrzK33PtfewOvbFsn6siFjksGwQAAPBfx86maFNkfJ7e43V7Nd09ZJ4spUor9fCvyjg1Ltc+NrtDmyLjFRmXomrB7rekEBNKAACA/5qzLVoWs0mSZM/K0KmpA3Vq6kDZszKy97Glp+jPiT11Zna4HHabzN6+spS6eUm0mE36Zqt7XktJoQQAAPivDUfispcHMnt6q9wzw2RNitX5jV9n75O4+gvZM9IU1GaYTGbLLR/bZndow9E4l2cuCjjlDQAAIOlihlXRiWk5XvO+M1T+DTvqwtYl8q3RSLbU80o7tFEBLfrJM/CuPH9GdEKaUjOsbveYRvf6agAAAPLpZEKqrnUrTtkmf1d6ZIQSvv9U9qxL8q54r0rXa5uvz3BIikpIVdidZZzKWtRwyhsAAEBS5nWWCTJZPBX09BBZk8/KkZmuoDZDZTKZXP45xRmFEgAAQJKXx/Vr0aUTv0uSHNZMWZNiC+xziiv3+4oAAADyoVKQn641d8yMO6Hzm+fJr05LeVWoqoSVE2W/lJqvzzD993PcDYUSAABAkp+3h0KuepKNw2ZVwg8TZLktSIEt+yuozVDZUpOUuG5avj4jJMjX7W7IkSiUAAAA2ZqFBmevQylJyVsWKPPscZV7eojM3r7yCq6sso88p9R9a5X+R0T2fuc3z9f5zfOVdmSLJOnigQ3Zr/3FYjapWY3gwvtiCpH7VWQAAIB86t4gRLN+i5IkZZyJVPJvC1X6oWdyPJfbv2EnpR3bqoSVE3Xni5Nl9rlNyZu+yXGc1L1rsv+77CPdJF1eh7JHw5CC/yIMYHI4HLf+sEoAAAA313P6Nm3+I155eJz3TVnMJjWuEqTZfRu47qBFCKe8AQAA/svhcKh22h7ZMjMkF87cPMwmjW1fx2XHK2oolAAAAJIyMjI0YMAAvTF4gB6wH5OcWGvyamPahqniVTf8uBMKJQAAKPFOnz6tZs2a6auvvtKMGTO0/JPXFd66hkuOPbJ1qLrWd89rJ//CNZQAAKBE27p1qzp06CCTyaSlS5eqQYP/Xec4PyJa76w4IKvdIVseLqq0mE3yMJs0pm2Y25dJiUIJAABKsBkzZuill17SQw89pCVLluiOO+7ItU9MYpreWLZPmyLjZTGbblgs/9retFo5jW1fx61Pc1+JQgkAAEqcrKwsDR8+XJMmTVK/fv00ceJEeXt73/A9x86maM62aG04GqfohDRdWaBMurxoebMawerRMETVgksXaP6ihkIJAABKlLi4OHXu3FlbtmzRpEmTNGDAgDwfIzXDqqiEVGVa7fLyMKtSkJ9bPgHnVlEoAQBAibFz5061b99emZmZWrx4sZo0aWJ0JLfAXd4AAKBEmDNnjpo0aaIKFSpox44dlEkXolACAAC3ZrVaFR4erh49eqhLly7auHGj7r77bqNjuZWSe7IfAAC4vYSEBHXr1k0bNmzQhAkTNHjwYJlcuGA5LqNQAgAAt7R37161a9dOFy5c0OrVq9W8eXOjI7ktTnkDAAC3s3jxYjVq1Ej+/v7asWMHZbKAUSgBAIDbsNlsGj16tDp37qy//e1v2rx5sypVqmR0LLfHKW8AAOAWzp8/rx49eujHH3/Uhx9+qJEjR3K9ZCGhUAIAgGLv0KFDateuneLi4vTjjz/qySefNDpSicIpbwAAUKytWLFCDRo0kIeHhyIiIiiTBqBQAgCAYslut+uf//ynnn32WbVo0UJbt25VtWrVjI5VInHKGwAAFDspKSnq1auXli1bpjFjxmj06NEym5mTGYVCCQAAipXIyEg9++yziomJ0fLly9W2bVujI5V4VHkAAFBs/PTTT6pfv76ysrK0bds2ymQRQaEEAABFnsPh0EcffaQ2bdqocePG2r59u2rVqmV0LPwXhRIAABRpqampeu655/T666/rH//4h1asWKGyZcsaHQtX4BpKAABQZEVFRaldu3aKjIzUokWL1KlTJ6Mj4RoolAAAoEhav369unTpIn9/f/3222+qU6eO0ZFwHZzyBgAARYrD4dBnn32m1q1bq27duoqIiKBMFnEUSgAAUGSkp6frhRde0NChQzV06FCtXLlSQUFBRsfCTXDKGwAAFAkxMTHq0KGD9u/fr2+++Ubdu3c3OhJuEYUSAAAY7tdff1XHjh3l7e2tzZs368EHHzQ6EvKAU94AAMBQX375pZo1a6aaNWtqx44dlMliiEIJAAAMkZGRoQEDBuill17SwIEDtXbtWgUHBxsdC/nAKW8AAFDoTp8+rU6dOmnHjh2aPn26+vTpY3QkOIFCCQAACtW2bdvUoUMHORwO/fLLL2rYsKHRkeAkTnkDAIBCM3PmTD366KO65557tHPnTsqkm6BQAgCAApeVlaXBgwerT58+ev7557VhwwbdcccdRseCi3DKGwAAFKhz586pc+fO2rx5s7744gsNGDBAJpPJ6FhwIQolAAAoML///rvat2+vS5cuaf369WratKnRkVAAOOUNAAAKxNy5c9WkSROVL19eO3bsoEy6MQolAABwKavVqpEjR6p79+7q1KmTNm3apIoVKxodCwWIU94AAMBlEhMT1a1bN61fv14TJkzQ4MGDuV6yBKBQAgAAl9i3b5/atWun5ORkrVq1Si1atDA6EgoJp7wBAIDTlixZokaNGum2225TREQEZbKEoVACAIB8s9vtevPNN9WpUye1adNGW7ZsUeXKlY2OhULGKW8AAJAvycnJ6t69u3788UeNGzdOr732GtdLllAUSgAAkGeHDx/Ws88+q7i4OP3444968sknjY4EA3HKGwAA5Ml3332nhx9+WB4eHtq+fTtlEhRKAABwa+x2u/75z3+qbdu2atGihbZu3arq1asbHQtFAKe8AQDATaWkpOiFF17Q0qVL9d577+nNN9+U2cxcCpdRKAEAwA1FRkaqXbt2io6O1vLly9W2bVujI6GI4Z8WAADgulatWqX69esrMzNT27Zto0zimiiUAAAgF4fDoX/96196+umn1ahRI23fvl21atUyOhaKKAolAADIIS0tTX//+9/12muvadSoUfruu+9UtmxZo2OhCOMaSgAAkC0qKkrt27fX0aNHtXDhQnXu3NnoSCgGKJQAAECStGHDBnXp0kWlS5fWb7/9pvvuu8/oSCgmOOUNAEAJ53A49Pnnn6tVq1a6//77FRERQZlEnlAoAQAowS5duqTevXtryJAhGjp0qH766ScFBQUZHQvFDKe8AQAoof7880916NBB+/bt0+zZs9WjRw+jI6GYolACAFACbd68WR07dpSXl5d+/fVXPfTQQ0ZHQjHGKW8AAIqZ1AyrDsQma1d0kg7EJis1w5qn90+ZMkXNmjVTaGioduzYQZmE05hQAgBQDBw7m6I526K14UicohPT5Lhim0lSSKCvmoUGq3uDEFWvUPqax8jMzNSrr76qqVOn6uWXX9ann34qT0/PQskP92ZyOByOm+8GAACMEJOYpjeW7dOmyHhZzCbZ7Nf/tv3X9qbVymls+zqqGOibve3MmTPq2LGjduzYocmTJ6tv376FER8lBIUSAIAian5EtN5ZcUBWu+OGRfJqFrNJHmaT3msbpm71Q7R9+3Z16NBBdrtdS5cuVcOGDQswNUoiCiUAAEXQpA3HNH71UaeP07xcquaOfl5169bVkiVLdOedd7ogHZATN+UAAJAPERERaty4sfz8/GQymbR79+7sbYMGDVKrVq3yfez5EdEuKZOStD7eT81eHK2ff/7ZZWVy1KhRatCggUuOBffAhBIAgDzKyspS9erV5ePjo+HDh8vX11dt2rRRQECATpw4odDQUK1atUrNmjXL87FfG/2uJi/6Semnjsiedl5lHnlOZZt2dyKtQ94eFq0d9liOayqvFhMToxkzZuiHH37QsWPHZLFYdO+99+rNN99Uy5Ytc+x75swZVapUSQsXLlTbtm2dyAZ3QaEEACCPDh8+rFq1amnatGl68cUXc2wbOnSoVq5cqSNHjuTr2CaTSRa/AHkGV9alE7+7oFBevqaycZUgze57/anipEmT9Nprr6ldu3Z65JFHZLVa9fXXX+v333/XjBkz1Lt37xz7d+3aVadPn9bGjRudygb3wLJBAADkUVxcnCSpbNmyOV7PysrSnDlzNHDgwHwd99jZFN01cLo8ylaQLS1Zf37uXJH8i83u0KbIeEXGpaha8LWXFGrWrJmio6NVrly57NcGDhyoBx54QG+//XauQtmlSxd17txZx48fV5UqVVySE8UX11ACAJAHL7zwgh577DFJUufOnWUymfT4449Lkn799VfFx8fnOkXcq1cv+fj46NChQzlef+KJJxQQEKDY2FhJ0pxt0fIOvN2lec9vmqOT4/6mzJN79M3W6OzX+/fvLy8vL+3Zs0eSFBYWlqNMSpK3t7eefvpp/fnnn0pJScmx7a+vcfny5S7Ni+KJQgkAQB4MGDBAb7zxhiRp8ODBmj17tkaPHi1J2rJli0wmk+rWrZvjPZ999pnKly+vXr16yWazSbr8tJrVq1dr4sSJ2TfLbDgSl6flgW5FmcZd5VWhiuJ++Exr9kZJklatWqVp06bp7bff1v3333/D9585c0a+vr7y9c15/WWZMmVUtWpVbd682aV5UTxRKAEAyINGjRpl38HdtGlT9ejRI/vnhw8fVmBgoPz9/XO8p2zZspo+fboiIiI0btw4nThxQuHh4WrXrp169OghSbqYYVV0YprL85osHgp6ZphsqYnau/hznTobr759+6pevXoaNWrUDd8bGRmppUuXqmPHjrJYLLm2V6lSRQcPHnR5ZhQ/FEoAAFwkISFBAQEB19zWunVrDRgwQGPGjFGHDh3k4+OjKVOmZG8/mZCqgrpL1qt8JZVt0l0X96zWU089qfj4eH311Vfy8Lj+rRRpaWnq3LmzSpUqpXHjxl1zn4CAAMXHxxdQahQnFEoAAFzoRounjB8/XoGBgdq9e7c+//xzBQcHZ2/LtNoLNJd/gw7yDK6sfbt26p133lHt2rWvu6/NZlO3bt108OBBLV68+LrrVzocDplMpoKKjGKEQgkAgIsEBQUpKSnputt37dqVfYf4vn37cmzz8ijYb8nW82dkTYq95mdfrV+/fvr+++81a9YsNW/e/Lr7JSUl5bqRByUThRIAABepWbOmkpKSlJycnGtbamqqevfurdq1a6t///766KOPFBERkb29UpCfCmrW53DYlfDDBJm9fBX+2ijNmzdPS5cuvea+I0eO1MyZM/Xpp5/queeeu+FxT5w4oVq1ahVEZBQzFEoAAFykUaNGcjgc2rlzZ65tr7/+uqKjo/XVV1/pk08+UaVKldSrVy9lZGRIkvy8PRRygyfZOOPC9m+VceqQanYZoQ8/eF+NGzfWSy+9lOv6x3/9618aP3683njjDQ0ZMuSGx0xOTtYff/yhxo0bF0hmFC8USgAAXKRJkyYKCgrS2rVrc7y+fv16TZ48WaNHj9aDDz4oPz8/zZw5U0eOHNFbb72VvV+z0GClHdig85vn60LEt5KkSzEHdH7zfJ3fPF/W5LjsfS+d3KuT457R+U1zbpgpKz5GyZu+0W33tVSndu1kNps1a9YsXbx4UYMGDcreb9myZXrttddUvXp11apVS998802OH2fPns1x3LVr18rhcOjZZ5/N7y8X3AhPygEAwEW8vLzUvXt3LVq0SGPHjpUkpaSkqE+fPqpbt272epXS5SWHhgwZoo8//lgdOnRQw4YN1b1BiMa9skoZMfuz98uI3quM6L2SJJ+7a8ujzOUbeexZlyRJltsCr5vHYbcp/odPZC7lr4Dm/dSjYYgkqXr16vrggw80ZMgQLVy4UF26dMle4PzYsWPq2bNnrmNt2LBBFSpUyP75okWL1KRJE1WtWjVfv1ZwLzzLGwAAFzp+/Lhq1qyplStXqkWLFnl+f8/p27T5j3jdbH3zpA0zlHpwo+4aME0mD88b7nsrz/LOizNnzqhy5cqaP38+E0pI4pQ3AAAuVaVKFfXt2/e6azfezMPm47JlZUg3mfdcOrlPZR7pdtMyKUkeZpPGtq+TrzzXMmHCBNWpU4cyiWxMKAEAKAKysrL02muvacKECWo18F0dLVvPZcf+sEMdda0f4rLjAVdjQgkAgMHi4uLUqlUrTZo0SZMmTdKqyW8rvHUNlxx7ZOtQyiQKHBNKAAAMtGPHDrVv316ZmZlavHixmjZtmr1tfkS03llxQFa7Q7abXVR5BYvZJA+zSWPahlEmUSiYUAIAYJBZs2apSZMmuuuuu/T777/nKJOS1K1+iNYOe0yNqwRJulwUb+Sv7Y2rBGntsMcokyg0TCgBAChkmZmZGjZsmCZPnqwXX3xRkyZNkre39w3fc+xsiuZsi9aGo3GKTkjTld+8TZJCgnzVrEawejQMUbXg0gWaH7gahRIAgEJ05swZderUSdu3b9ekSZPUv3//PB8jNcOqqIRUZVrt8vIwq1KQn/y8WVoaxqFQAgBQSLZu3aqOHTvK4XBoyZIlatSokdGRAJfgGkoAAArB1KlT9eijj6py5crauXMnZRJuhUIJAEABysjIUP/+/TVgwAD169dP69ev1x133GF0LMCluOACAIACcurUKXXs2FG7d+/WjBkz1Lt3b6MjAQWCQgkAQAHYtGmTOnfuLE9PT23atEn169c3OhJQYDjlDQCACzkcDv373/9W8+bNVbNmTe3cuZMyCbdHoQQAwEXS09PVu3dvvfLKK3r55Ze1Zs0aBQcHGx0LKHCc8gYAwAWio6PVoUMHHThwQLNnz1aPHj2MjgQUGgolAABO2rBhg7p06SI/Pz9t2bJFdevWNToSUKg45Q0AQD45HA59+umnatWqle6//37t2LGDMokSiUIJAEA+pKWlqUePHho+fLiGDx+un376SeXKlTM6FmAITnkDAJBHJ06cUPv27XXs2DHNnz9fXbt2NToSYCgKJQAAebBmzRp169ZNZcuW1datW1WnTh2jIwGG45Q3AAC3wOFw6KOPPtKTTz6phx9+WBEREZRJ4L8olAAA3MTFixfVtWtXvf766xo1apS+//57BQYGGh0LKDI45Q0AwA1ERkaqffv2ioqK0pIlS9ShQwejIwFFDhNKAACu48cff1T9+vWVmZmpbdu2USaB6yjxhTI1w6oDscnaFZ2kA7HJSs2wGh0JAGAwu92u//u//9Mzzzyjpk2bavv27apdu7bRsYAiq0Se8j52NkVztkVrw5E4RSemyXHFNpOkkEBfNQsNVvcGIapeobRRMQEABrhw4YJ69eqlb7/9Vu+++67eeustmc0lfv4C3JDJ4XA4br6be4hJTNMby/ZpU2S8LGaTbPbrf+l/bW9arZzGtq+jioG+hZgUAGCEI0eOqF27doqNjdU333yjv/3tb0ZHAoqFElMo50dE650VB2S1O25YJK9mMZvkYTbpvbZh6lY/pAATAgCMtHz5cvXs2VN33323li1bptDQUKMjAcVGiZjhT9pwTKOW7lOG1Z6nMilJNrtDGVa7Ri3dp0kbjhVQQgCAUex2u95++221a9dOrVq10rZt2yiTQB65faGcHxGt8auPuuRY41cf1YKIaJccCwBgvPPnz6tt27b6v//7P73//vtavHixSpfm2nkgr9z6lHdMYppafvqLMqz2m+57KXq/Lmxfqsyzx2VLS5bZx09ewVVU5pFu8rn7f3f2eXuYtXbYY1xTCQDF3IEDB9S+fXudO3dOc+fO1VNPPWV0JKDYcusJ5RvL9sl6i6e4s5JOSSazStd9SoGtX5L/wx1kS03S2TmjlH58Z/Z+VrtDbyzbV1CRAQCFYPHixWrQoIF8fHy0Y8cOyiTgJLedUB47m6JWEzY6dQx71iWd+vJFeQVXUYWuY3JsWzvsUVUL5rQIABQnNptNb775psaNG6euXbtq+vTp8vPzMzoWUOy57YRyzrZoWcwm2bMydGrqQJ2aOlD2rIzs7bb0FP05safOzA6Xw2675jHMnj6y+JaRPSM1x+sWs0nfbOVaSgAoThITE9WmTRt99NFH+te//qV58+ZRJgEXcdtCueFInGx2h8ye3ir3zDBZk2J1fuPX2dsTV38he0aagtoMk8lsyX7dnpEmW1qyshJilPTLV8o6d1I+99yf49g2u0MbjsYV2tcCAHDOnj17VK9ePUVERGjVqlUKDw+XyWQyOhbgNtzySTkXM6yKTkzL/rn3naHyb9hRF7YukW+NRrKlnlfaoY0KaNFPnoF35XjvuW/H6dKJ3y//xOKh2x54UmUf6ZbrM6IT0pSaYZWft1v+EgKA25g3b5769u2r0NBQrV+/XpUqVTI6EuB23LINnUxI1dUXhpZt8nelR0Yo4ftPZc+6JO+K96p0vba53hvw+AuyPdxetgvndHH/ejlsVjnsNl3971iHpKiEVIXdWaagvgwAgBOsVqtGjRqljz/+WN27d9fUqVPl68sKHUBBcMtT3pnXWCbIZPFU0NNDZE0+K0dmuoLaDL3m6Q6vClVUqnJd3XZ/a1Xo9k9lnj6qhB8+veXPAQAYLz4+Xk888YQmTJigCRMmaPbs2ZRJoAC5ZaH08rj2l/XXqWyHNVPWpNibHsdk8VSp6g2UduS3HDf0/OX0nzGyWq3OhQUAuNTvv/+uhx56SPv27dPatWs1ZMgQrpcECphbFspKQX65TlFnxp3Q+c3z5FenpbwqVFXCyomyX0q95vuv5MjKlOSQIzM95+sOh555vIH8/Px0//33q3v37ho7dqyWL1+uyMhI2WzXvnMcAFBwZs+erUceeUTBwcHauXOnHn/8caMjASWCW15D6eftoZBAX5387405DptVCT9MkOW2IAW27C9r8lmd/mqYEtdNU7k2QyVJttTzsviVzXEc+6WLSju6RRb/8rm2VSzro5k//aD9+/frwIEDOnDggH788UedP39eklSqVCnVqlVLYWFhuvfeexUWFqawsDCFhITIbHbLHg8AhsnKylJ4eLg+//xz9e7dW5MnT5aPj4/RsYASwy0LpSQ1Cw3W7G0nZbM7lLxlgTLPHleF596X2dtXXsGVVfaR53R+42z51XxEparWV9zCd2QpXU7ed4bK7FdGtuRzurhvrWwXE1Xu2ddyHNtiNqll7TvUvHmYmjdvnv26w+FQbGxsdsH8q2wuW7ZMFy9elCT5+flll8sry+Zdd93FKRkAyIezZ8+qS5cu2rJliyZPnqyBAwfy9ylQyNz+STkZZyJ15usRKl33aQW2GpC93WG36czscNlSEnTni5OVeuBnpR7aqKyEP2XPSJXZ57bLyw016CCfivfmOn5enpTjcDgUExOTY5p54MABHTx4UGlpl6eoZcqUUe3atXNMM++9915VqFCBvxgB4Dq2b9+uDh06yGazadGiRWrSpInRkYASyW0LpST1nL5NW44nyHaLz/O+FRazSY2rBGl23wZOH8tutysqKirHNPPAgQM6dOiQMjIu3wQUGBiY67R5WFiYypcv7/TnA0BxNn36dA0aNEgPPviglixZojvvvNPoSECJ5daFMiYxTS0//UUZLlzex9vDrLXDHlPFwIJbfsJms+mPP/7Ider8yJEjysrKkiQFBwfnOm0eFhamgICAAssFAEVBZmamhgwZoi+//FIDBgzQZ599Jm9vb6NjASWaWxdKSZofEa1RS/e57HgfdqijrvVDXHa8vMjKylJkZGSOaeb+/ft17Nix7LvK77jjjlynzWvXri1/f39DMgOAK8XGxqpz587asWOH/v3vf+vFF180OhIAlYBCKUmTNhzT+NVHnT7OyNaherlZNRckcq2MjAwdPXo016nzyMhI/fW/t2LFirmmmbVr15afn5/B6QHg1mzZskUdO3aUxWLRkiVL1KCB85ceAXCNElEopcuTyndWHJDV7sjTNZUWs0keZpPGtA0zbDKZX+np6Tp8+HCuU+cnTpzI3qdy5cq5Tp3XrFlTpUqVMjA5APyPw+HQlClTNHjwYDVo0ECLFi3S7bffbnQsAFcoMYVSunxN5RvL9mlTZLwsZtMNi+Vf25tWK6ex7esU6DWThe3ixYs6dOhQjjvO9+/fr5iYGEmS2WxW1apVcy1vVKNGDa5TkpSaYVVUQqoyrXZ5eZhVKchPft5uuwIXYKhLly7p5Zdf1owZM/TKK6/o448/lpeXl9GxAFylRBXKvxw7m6I526K14WicohPSdOUvgElSSJCvmtUIVo+GIbe8NJA7uHDhgg4ePJhreaPY2MuPqbRYLKpevXquO86rV68uT09Pg9MXrOzfM0fiFJ14jd8zgb5qFhqs7g1CVL1Cyfk9AxSkmJgYdezYUXv37tWUKVPUq1cvoyMBuI4SWSivxLTp5pKSknKdNj9w4IDi4uIkSZ6engoNDc11jWbVqlVlsVgMTu8cptqAMTZu3KjOnTvL29tby5Yt00MPPWR0JAA3UOILJfLv3LlzOSaZfxXOxMRESZK3t3f24yevLJuVKlUqFo+fdPa62/fahqlbMbvuFjCaw+HQxIkTNWLECDVp0kQLFy5k3V2gGKBQwqUcDofOnj2ba5p54MABJScnS5J8fX1Vq1atXMsbVaxYscg8FchVKwOEt66hV5pVd0EiwP2lp6drwIABmj17toYPH64PP/xQHh6cMQKKAwolCoXD4dCpU6dynTo/ePBg9nPOS5cunevxk2FhYbrzzjtzFc2IiAgNGTJEe/bsUVpamnbt2qUHHnhAkjRo0CAdO3ZMa9asyVfWor526ahRo7RhwwZt27bNZccEjBYVFaUOHTro8OHD+s9//qO///3vRkcCkAcUShjKbrfnes75/v37dejQIaWnp0uSypYtm2OSGRoaqr59+8rPz0/Dhw+Xr6+v2rRpo4CAAJ04cUKhoaFatWqVmjVrlqcshw8f1oTJUzRr4QplJp2W2dNHXrdXVZkm3eV9R/6njLfydKX09HS98sor2rZtm2JiYmSz2VS1alX16dNHgwYNynHT05kzZ1SpUiUtXLhQbdu2zXcuoKhYt26dunbtqtKlS+vbb7/V/fffb3QkAHlEoUSRZLPZrvuc88zMTEmSn5+f6tWrl6Nszp07Vxs2bNCRI0fy/Jnh4eGa+MVUeVdvJM87asiekaqLu36SNfmsgruOUalKD+Tra7mV578nJibq6aef1qOPPpp9jemWLVv0zTffqFu3bpo7d26O/bt27arTp09r48aN+coEFAUOh0OffPKJXnvtNbVo0ULz5s1TUFCQ0bEA5AOFEsWK1WrV/Pnz1bNnT3Xt2lUOhyP7OedWq1XS5aLZoEGDXKfOy5Yte8NjL1u9UUNWxcns9b9F3W3pFxQ77SV5Bt6l23t85FT2tcMezfMyVK+++qomTZqk06dP51jIecmSJercubMiIyNVpUoVp3IBRkhNTdWLL76o+fPn6/XXX9f7779f7FeFAEqyon+rLXCFF198UT179pQkLViwQAsXLlS5cuWUmpqqGTNmSJI6duyogIAArVq1SoMGDVKTJk0UEBCgChUq6IknntDw4cM1Y8YMNWzYUAEBAdnrbO65FCRPn5ynpS2l/OVzd5iy4mPylff8pjk6Oe5vyjy5R99sjc5+vX///vLy8tKePXtu+P5KlSpdPs758zleb9mypSRp+fLl+coFGOn48eNq1KiRvvvuOy1cuFDjxo2jTALFHLfPoVgZMGCA7rrrLo0dO1aDBw9W/fr1VaFCBXl5eSk2NlYmk0kTJ06Uv7+/pMvPOY+IiNCzzz4ri8UiHx8ffffdd/r000+zj9m4cWOFhYUpqk5v2Uy5HzlpS02S2dc/X3nLNO6q9MjtivvhM62pfb/ebRumVatWadq0afrnP/+Z61qxzMxMXbhwQenp6dqxY4fGjx+ve+65R9Wq5XyGfJkyZVS1alVt3rxZw4YNy1c2wAirVq3Sc889p8DAQG3dulX33nuv0ZEAuAATShQrjRo1UqtWrSRJTZs2VY8ePbJ/fvjwYQUGBmaXSenyWphNmjTRvHnzdPr0aT388MNavXq1brvtNj3++OP6+uuv1a1bN9ktXkqVT67PuxSzXxmnDsuvZtN85TVZPBT0zDDZUhO1d/HnOnU2Xn379lW9evU0atSoXPsvXbpU5cuXV0hIiDp06KC7775b33333TWXTqlSpYoOHjyYr1xAYXM4HPrggw/01FNPqVGjRoqIiKBMAm6ECSXcRkJCggICAq65rXXr1howYIDGjBmjxYsXy8fHRwsWLFBwcLAk6UBsstpM/DXHe2yp5xW/Yrw8ylaQf8OO+c7lVb6SyjbprvO/fKWnnnpS8fHxWr169TVLYrNmzbRmzRqdP39e69at0549e5SamnrN4wYEBGjXrl35zgUUlpSUFPXu3VtLlizRW2+9pXfffbdYPNwAwK2jUMKt3Oges/Hjx2v58uXavXu35s6dm10mJSnTas+xrz3zkuIWvyd7Zrpu7/Fhjht18sO/QQelHtqofbt2auzYsapdu/Y196tQoYIqVKggSerUqZPGjh2rVq1a6dixYzluypEuf61FZSF44HqOHTumdu3aKSYmRsuWLVO7du2MjgSgAPBPRLiNoKAgJSUlXXf7rl27sp8/vm9fzoXLvTz+90fBYcvSuWXvKzMuSsEd35RX+UpOZ7OePyNr0uWbf3788UcdPHgwe/mjG+nUqZMuXrx4zZtvkpKSVK5cOaezAQXl+++/V/369WWz2bRt2zbKJODGKJRwGzVr1lRSUlL2Ix6vlJqaqt69e6t27drq37+/PvroI0VERGRvrxTkJ5Mkh8Ou+O8/0aWoPSrfdqR8Quo4ncvhsCvhhwkye/nKu5Svfv31V4WFhcnX11ehoaF69tln9frrr2vmzJn67bffcpTivxZ3v9bXdOLECdWqVcvpfICr2e12jRkzRn/729/0+OOPa/v27fxeBdwcp7zhNho1aiSHw6GdO3eqefPmOba9/vrrio6O1tatWxUaGqp169apV69e2rVrl7y9veXn7aGQQF/9Pv9jpR3apMAnX5FvaGOX5Lqw/VtlnDqksN7/pz3TRqlRo0aKjIzU6NGj9eeff+rQoUOaO3eu/vzzz+z3VKhQQTVr1lR8fLyky9+go6KiFBISIrPZrOTkZP3xxx966aWXXJIRcJXk5GQ9//zz+u677zRmzBiNHj2a6yWBEoBCCbfRpEkTBQUFae3atTkK5fr16zV58mS98847evDBByVJM2fO1OOPP6633npLH310ecHyUkdW6eLvP8j7rpoyeXjr4v4NOY7vW6ORzF6X7wS/dHKvzs57Q2UeeU5lm3a/bqas+Bglb/pGt93XUp3atZPFYtGcOXP0wAMPaNu2bVq4cKEkacKECfriiy/UuHFjeXh46OTJk9q7d6/Onj0rs9msf/zjH/rHP/6hUqVKKTQ0VH5+fnI4HLJYLNqzZ49q1KihUqWcu84TcNahQ4fUvn17nTlzRt99953atGljdCQAhYRCCbfh5eWl7t27a9GiRRo7dqyky3eX9unTR3Xr1tXo0aOz923atKmGDBmijz/+WB06dFDDhg3lm3J58fKMU4eVcepwruP7DJyeXSjtWZckSZbbAq+bx2G3Kf6HT2Qu5a+A5v3Uo2GIJKl69er64IMPNGTIEC1cuFBdunRRkyZNtGXLFq1fv15nz56Vh4eHQkND9frrr2vQoEGKjY3V4cOHdfjwYR06dEjLli2Th4eHhg4dKkkymUyqVKmSatasmetH+fLluXkHBW7ZsmV6/vnnFRISooiICFWvXt3oSAAKEY9ehFs5fvy4atasqZUrV6pFixZ5fn/P6du05XiCbPYb/7FI2jBDqQc36q4B02Ty8LzhvrfyLO+8OHPmjCpXrqz58+erSZMmOnLkSHbR/Kt0Hj9+XHb75TvXAwMDr1k0K1eufM2li4C8sNlseuedd/T++++rY8eOmjlzpkqXztsjRgEUfxRKuJ2XXnpJkZGRWrNmTZ7fG5OYppaf/qKMq5YRutrpWcN02wNPqPQDT970mN4eZq0d9pgqBvredN9bMWrUKK1fv17bt2+/7j4ZGRmKjIzMVTQPHz6cva6ll5eXqlevnqtohoaGUghwS5KSktS9e3f99NNPGjt2rF5//XWm4UAJRaEErjI/Ilqjlu67+Y636MMOddS1fojLjucMh8OhU6dO5Th9/td///VMc0m6++67rznVvPPOOykMkCTt379f7dq1U2JioubNm6cnnnjC6EgADEShBK5h0oZjGr/6qNPHGdk6VC83q3bzHYuACxcuXPP0+bFjx2S1WiVJpUuXvmbRrFatmry8vAz+ClBYFi5cqN69e6tatWpatmyZqlSpYnQkAAajUALXMT8iWu+sOCCr3XHTayqvZDGb5GE2aUzbsCIzmXRGVlaWTpw4kWuqeejQoez1MS0Wi6pWrXrNsnm9x2Gi+LFarRo9erQ++ugjPffcc5o2bZr8/PyMjgWgCKBQAjcQk5imN5bt06bIeFnMphsWy7+2N61WTmPb13HZNZNFlcPhUFxc3DVPn588eTJ7v7/W1Lz6x19raqJ4SEhIULdu3bR+/Xr961//0rBhw7j8AUA2CiVwC46dTdGcbdHacDRO0QlpuvIPjUlSSJCvmtUIVo+GIaoWzA0taWlpOnr0aK6ieeTIEWVkZEhS9pqaVxdN1tQsenbv3q327dvr4sWLWrBgQa4HBwAAhRLIo9QMq6ISUpVptcvLw6xKQX7y82b5nVths9kUHR19zanmuXPnJLGmZlEzZ84c9evXT7Vq1dLSpUt1zz33GB0JQBFEoQRQJCQkJLCmpos5848fq9WqkSNHasKECXr++ef15ZdfMjkGcF0USgBFGmtq5k325RlH4hSdeI3LMwJ91Sw0WN0bhKh6hWv/2sTFxalr16769ddf9cknn+iVV15hMgzghiiUAIol1tTMyVU3kO3YsUMdOnRQRkaGFi1apEcffbQw4gMo5iiUANxOSVtT09klrt5rG6Zu9UM0a9YsDRw4UPfff7+WLFmiu+++uwBTA3AnFEoAJYY7rqnpqkX4QzOOavWnw9W3b19NmjRJPj4+LkgHoKSgUAIo8Yrrmpqufkzo00GJ+veIHm51KQCAwkGhBIAbKKprasYkpqnlp78ow2q/6b7pUbuVeuBnZfx5ULaUBFn8ysrnnvtV5tEe8rgtMHs/bw+z1g57zO0X5QfgehRKAMgHo9fU7Dl9m7YcT7ilayZPzxoqe3qKfGs2kUfgnbKeP6OUnd/L5OmtO3tPlOW2y6fyLWaTGlcJ0uy+DfKdC0DJRKEEABcr6DU1j51NUasJG285z6Xo/fKuWFsmkznHa2fnjpJ/464KeLRnjv3XDnuUJz4ByBMKJQAUEletqfnuigOave2ksjIu6fTMIZKkO3p/JrOntyTJlp6i0/8ZJI+yFVSh+4cymS3XzBMz4Tn5hNRR+Q5vZL9mMZvUs8E9erdtWEH+UgBwMzxaAgAKibe3t8LCwhQWlrOsXW9NzZkzZ15zTc1TD/aTzewns6e3yj0zTGdmj9T5jV8rsEU/SVLi6i9kz0hTUJth1y2T9sx02bPSZfb1z/G6ze7QhqNxelcUSgC3jkIJAAYzmUy6++67dffdd6tly5Y5tl29puaBI38ozfS/m2a87wyVf8OOurB1iXxrNJIt9bzSDm1UQIt+8gy867qfeSFiuWSzyrdm01zbohPSlJph5Rn1AG4Zp7wBoBg5EJusNhN/zfGaw5al07OGyZGZLnvWJXkGVVSFv39w3Zt+LkXv19n5o+Vbo7HKt3v9mvv88GoThd1ZxuX5Abinwl84DQCQb5nXWCbIZPFU0NNDZE0+K0dmuoLaDL1umcxKiNG5pe/Ls9w9Cnrq1Tx9DgBcD4USAIoRL49r/7V96cTvkiSHNVPWpNhr7mO9cE5nF7wts7evgru8K7P39debvN7nAMC18DcGABQjlYL8dPXsMTPuhM5vnie/Oi3lVaGqElZOlP1Sao59bOkXdHbBW3JYsxTcdUyOBc2vZvrv5wDAraJQAkAx4uftoZArnmTjsFmV8MMEWW4LUmDL/gpqM1S21CQlrpuWvY8985LiFr4rW0qCgru8e8ObdSQpJMiXG3IA5AmFEgCKmWahwbKYL88pk7csUObZ4yr39BCZvX3lFVxZZR95Tqn71ir9jwhJUvx345V5+qh8azZRVnyMLu7fkP0j7ehvOY5tMZvUrEZwoX9NAIo3/gkKAMVM9wYhmvVblDLORCr5t4Uq/dAz8rnnvuzt/g07Ke3YViWsnKg7X5yszLPHJUmpe9code+aHMey+AfLt0aj7J/b7A71aBhSOF8IALfBskEAUAzl5Vnet4pneQPIL055A0AxNLZ9HXmYr700UH55mE0a276OS48JoGSgUAJAMVQx0Ffvufh522Pahqli4PWXEgKA66FQAkAx1a1+iMJb13DJsUa2DlXX+lw7CSB/uIYSAIq5+RHRemfFAVntjjxdU2kxm+RhNmlM2zDKJACnUCgBwA3EJKbpjWX7tCkyXhaz6YbF0uSwy2Ey65EqgRrX8X5OcwNwGqe8AcANVAz01ey+DbRm6KPq2eAe3RPkm+uJOiZJ9wT56m+1yurUtIFq7XWUMgnAJZhQAoCbSs2wKiohVZlWu7w8zKoU5Jf9BJxnnnlG0dHR2rNnj0wm194tDqDkoVACQAn0888/q1mzZlq1apVat25tdBwAxRyFEgBKIIfDoXr16ikoKEirV682Og6AYo5rKAGgBDKZTAoPD9eaNWu0d+9eo+MAKOaYUAJACZWVlaVq1arp8ccf11dffWV0HADFGBNKACihPD09NWTIEM2dO1enTp0yOg6AYoxCCQAl2IsvvihfX19NnDjR6CgAijEKJQCUYP7+/urfv7++/PJLpaSkGB0HQDFFoQSAEm7w4MFKTU3VjBkzjI4CoJjiphwAgHr27KlNmzYpMjJSHh4eRscBUMwwoQQAaMSIETp58qSWLl1qdBQAxRATSgCAJKlly5a6cOGCtm3bxuMYAeQJE0oAgKTLU8qIiAht2rTJ6CgAihkmlAAASZcfx3jvvfeqWrVqWr58udFxABQjTCgBAJIuP45xxIgRWrFihY4cOWJ0HADFCBNKAEC2jIwM3XPPPWrXrp2+/PJLo+MAKCaYUAIAsnl7e2vw4MH66quvFBcXZ3QcAMUEhRIAkMPAgQNlNps1efJko6MAKCYolACAHAIDA9WnTx/9+9//Vnp6utFxABQDFEoAQC5Dhw5VQkKCvv76a6OjACgGuCkHAHBNnTp10r59+3To0CGZzcwfAFwff0MAAK4pPDxcR48e1ffff290FABFHBNKAMB1PfLII/Lw8NAvv/xidBQARRgTSgDAdYWHh2vjxo3avn270VEAFGFMKAEA12Wz2RQaGqqHHnpICxYsMDoOgCKKCSUA4LosFouGDx+uxYsXKyoqyug4AIooCiUA4IZeeOEFlS1bVhMmTDA6CoAiikIJALghX19fDRo0SP/5z3+UlJRkdBwARRCFEgBwU6+88oqysrI0depUo6MAKIK4KQcAcEtefPFFrVy5UidOnJCXl5fRcQAUIUwoAQC3ZPjw4YqNjdX8+fONjgKgiGFCCQC4ZW3atFFMTIz27Nkjk8lkdBwARQQTSgDALQsPD9e+ffu0du1ao6MAKEKYUAIAbpnD4VC9evVUrlw5rVq1yug4AIoIJpQAgFtmMpk0YsQIrV69Wnv37jU6DoAiggklACBPsrKyVLVqVTVv3lyzZs0yOg6AIoAJJQAgTzw9PTVkyBDNnTtXsbGxRscBUARQKAEAedavXz+VKlVKEydONDoKgCKAQgkAyDN/f3/1799fX375pS5evGh0HAAGo1ACAPJl8ODBunjxombMmGF0FAAG46YcAEC+9ejRQ5s3b9axY8fk4eFhdBwABmFCCQDItxEjRigqKkpLly41OgoAAzGhBAA4pUWLFkpJSdG2bdt4HCNQQjGhBAA4ZcSIEYqIiNCvv/5qdBQABmFCCQBwit1uV506dVS9enV9++23RscBYAAmlAAAp5jNZo0YMUIrVqzQkSNHjI4DwABMKAEATsvIyNA999yjdu3a6csvvzQ6DoBCxoQSAOA0b29vvfrqq/rqq6907tw5o+MAKGQUSgCASwwcOFAmk0mTJ082OgqAQkahBAC4RFBQkPr06aN///vfSk9PNzoOgEJEoQQAuMywYcMUHx+v2bNnGx0FQCHiphwAgEt17NhRBw4c0MGDB2U2M7cASgL+pAMAXCo8PFxHjhzRDz/8YHQUAIWECSUAwOUaN24sT09P/fLLL0ZHAVAImFACAFwuPDxcGzduVEREhNFRABQCJpQAAJez2WwKDQ1VvXr1NH/+fKPjAChgTCgBAC5nsVg0bNgwLV68WFFRUUbHAVDAKJQAgALxwgsvqEyZMvrss8+MjgKggFEoAQAFws/PTy+99JL+85//6Pz580bHAVCAKJQAgALzyiuvKDMzU1OnTjU6CoACxE05AIAC1bdvX/300086ceKEvLy8jI4DoAAwoQQAFKjhw4crNjZWCxYsMDoKgALChBIAUODatGmjP//8U7t375bJZDI6DgAXY0IJAChwI0aM0N69e7Vu3TqjowAoAEwoAQAFzuFw6KGHHlJwcLB++ukno+MAcDEmlACAAmcymTRixAitWrVK+/btMzoOABdjQgkAKBRZWVmqUqWKWrZsqZkzZxodB4ALMaEEABQKT09PDRkyRHPmzFFsbKzRcQC4EIUSAFBo+vXrJx8fH02aNMnoKABciEIJACg0ZcqUUf/+/fXFF1/o4sWLRscB4CIUSgBAoRo8eLBSUlI0Y8YMo6MAcBFuygEAFLru3bvrt99+09GjR+Xh4WF0HABOYkIJACh0I0aM0IkTJ7Rs2TKjowBwASaUAABDNG/eXKmpqdq6dSuPYwSKOSaUAABDhIeHa/v27dq8ebPRUQA4iQklAMAQdrtd9957r0JDQzn1DRRzTCgBAIYwm80aMWKEli9frqNHjxodB4ATKJQAAMN0795d5cuX16effmp0FABOoFACAAzj4+OjV199VbNmzdK5c+eMjgMgnyiUAABDvfTSSzKZTPriiy+MjgIgn7gpBwBguJdfflmLFi1SdHS0fHx8jI4DII+YUAIADDds2DDFx8dr9uzZRkcBkA9MKAEARUKHDh108OBBHTx4UGYz8w6gOOFPLACgSAgPD9eRI0f0448/Gh0FQB4xoQQAFBmNGjWSt7e3fv75Z6OjAMgDJpQAgCIjPDxcv/zyi3bs2GF0FAB5wIQSAFBk2Gw21ahRQw8//LDmzZtndBwAt4gJJQCgyLBYLBo2bJgWLVqkkydPGh0HwC2iUAIAipTevXvL399fn332mdFRANwiCiUAoEjx8/PTSy+9pGnTpun8+fNGxwFwCyiUAIAi55VXXlFmZqamTZtmdBQAt4CbcgAARVKfPn20evVqHT9+XF5eXkbHAXADTCgBAEXSiBEjdOrUKS1cuNDoKABuggklAKDIevrppxUbG6tdu3bJZDIZHQfAdTChBAAUWSNGjNCePXu0bt06o6MAuAEmlACAIsvhcOjBBx/U7bffrpUrVxodB8B1MKEEABRZJpNJI0aM0E8//aT9+/cbHQfAdTChBAAUaVlZWapcubJat26tGTNmSJJSM6yKSkhVptUuLw+zKgX5yc/bw+CkQMnFnz4AQJHm6empoUOH6u3xkxU0f7u2x6QqOjFNV05DTJJCAn3VLDRY3RuEqHqF0kbFBUokJpQAgCItJjFNIxft0tao8zLJIYeuf7e3xWySze5Q02rlNLZ9HVUM9C3EpEDJRaEEABRZ8yOi9c6KA7LaHbLZb/3blcVskofZpPfahqlb/ZACTAhAolACAIqoSRuOafzqo04fJ7x1Db3SrLoLEgG4Hu7yBgDkS0REhBo3biw/Pz+ZTCbt3r07e9ugQYPUqlWrfB97fkS0S8qkJI1ffVQLIqJdcqy/dOvWTV26dHHpMYHijAklACDPsrKyVL16dfn4+Gj48OHy9fVVmzZtFBAQoBMnTig0NFSrVq1Ss2bN8nTcw4cPa8LkKZq1cIUyk07L7Okjr9urqkyT7vK+I/9TRm8Ps9YOeyxP11T++uuvatq0qSTp3LlzKleuXPa2Xbt2qV69evr99991//335zsX4C4olACAPDt8+LBq1aqladOm6cUXX8yxbejQoVq5cqWOHDmS5+OGh4dr4hdT5V29kTzvqCF7Rqou7vpJ1uSzCu46RqUqPZCvvBazSY2rBGl23wa3tL/dbtdDDz2kY8eOKTU1NVehlKQGDRooNDRUX3/9db4yAe6EU94AgDyLi4uTJJUtWzbH61lZWZozZ06+Twc/0rqtKgycocCnBqv0A0+qTIOOur3XxzKXKq3kX+fmO6/N7tCmyHhFxqXc0v5Tp05VTExMrrJ8pS5dumjp0qW6ePFivnMB7oJCCQDIkxdeeEGPPfaYJKlz584ymUx6/PHHJV0+TRwfH6+WLVvmeE+vXr3k4+OjQ4cO5Xj9iSeeUEBAgGJjYyVJey4FydMn52lpSyl/+dwdpqz4mHzlPb9pjk6O+5syT+7RN1v/dy1l//795eXlpT179uTYPzExUW+++abGjBmTqzBfqVWrVkpNTdWaNWvylQtwJxRKAECeDBgwQG+88YYkafDgwZo9e7ZGjx4tSdqyZYtMJpPq1q2b4z2fffaZypcvr169eslms0mSpkyZotWrV2vixIm68847JUkbjsRdc3kgW2qSzL7++cpbpnFXeVWoorgfPtOavVGSpFWrVmnatGl6++23c10D+dZbb+n222/XgAEDbnjc2rVrq1SpUtq8eXO+cgHuhEIJAMiTRo0aZd/B3bRpU/Xo0SP754cPH1ZgYKD8/XOWv7Jly2r69OmKiIjQuHHjdOLECYWHh6tdu3bq0aOHJOlihlXRiWm5Pu9SzH5lnDosv5pN85XXZPFQ0DPDZEtN1N7Fn+vU2Xj17dtX9erV06hRo3Lsu3fvXk2ZMkWffPKJLBbLDY/r4eGhihUr6uDBg/nKBbgTHr0IAHCZhIQEBQQEXHNb69atNWDAAI0ZM0aLFy+Wj4+PpkyZkr39ZEKqrp5N2lLPK37FeHmUrSD/hh3zncurfCWVbdJd53/5Sk899aTi4+O1evVqeXjk/DY4ePBgPfXUU2rduvUtHTcgIEDx8fH5zgW4CwolAMClbrR4yPjx47V8+XLt3r1bc+fOVXBwcPa2TKs9x772zEuKW/ye7Jnpur3HhzJ7lXIql3+DDko9tFH7du3U2LFjVbt27RzbFyxYoC1btmj//v23fEyHwyGT6fqPggRKCk55AwBcJigoSElJSdfdvmvXruw7xPft25djm5fH/74lOWxZOrfsfWXGRSm445vyKl/J6WzW82dkTYq95mdL0siRI9W5c2d5eXkpKipKUVFROn/+vCQpJiYm+8ahKyUlJeVaTggoiSiUAACXqVmzppKSkpScnJxrW2pqqnr37q3atWurf//++uijjxQREZG9vVKQn0ySHA674r//RJei9qh825HyCanjdC6Hw66EHybI7OWr8NdGad68eVq6dGmOfWJiYjR37lxVrlw5+8dnn30mSXrwwQf19NNP59jfarUqJiZGtWrVcjofUNxxyhsA4DKNGjWSw+HQzp071bx58xzbXn/9dUVHR2vr1q0KDQ3VunXr1KtXL+3atUve3t7y8/ZQSKCvfp//sdIObVLgk6/IN7SxS3Jd2P6tMk4dUljv/9OHH/xDW37dqJdeekmPPvpo9oRx2bJlud43f/58LViwQF9//bXuvvvuHNsOHjyoS5cuqXFj12QEijMKJQDAZZo0aaKgoCCtXbs2R6Fcv369Jk+erHfeeUcPPvigJGnmzJl6/PHH9dZbb+mjjz6SJJU6skoXf/9B3nfVlMnDWxf3b8hxfN8ajWT28pEkXTq5V2fnvaEyjzynsk27XzdTVnyMkjd9o9vua6lO7drJbDZr1qxZeuCBBzRo0CAtXLhQktSuXbtc7/3r+eRPPfVUrlPba9aska+vr1PPLAfcBae8AQAu4+Xlpe7du2vRokXZr6WkpKhPnz6qW7du9nqV0uUlh4YMGaKPP/5YW7dulST5plxevDzj1GElfP9xrh/2tP+dSrdnXZIkWW4LvG4eh92m+B8+kbmUvwKa91OPhiGSpOrVq+uDDz7QokWLsgtlXi1atEgdOnRQ6dKl8/V+wJ3wLG8AgEsdP35cNWvW1MqVK9WiRYs8v7/n9G3acjzhmgucXylpwwylHtyouwZMk8nD84b75vVZ3jeze/duPfjgg/r999/1wAMPuOSYQHHGhBIA4FJVqlRR3759NW7cuHy9f2z7OvIw33wpnksn96nMI91uWiYlycNs0tj2zt/c85dx48apU6dOlEngv5hQAgCKnPkR0Rq1NPfSPvn1YYc66lo/xGXHA5ATE0oAQJHTrX6IwlvXcMmxRrYOpUwCBYwJJQCgyJofEa13VhyQ1e646TWVV7KYJA+LWWPahlEmgUJAoQQAFGkxiWl6Y9k+bYqMl8VsumGxtJgkm0O6w5SsheHtVDHQtxCTAiUXhRIAUCwcO5uiOduiteFonKIT0nTlNy+TpJAgXzWrEayT677Rd3OmKTo6Wr6+FEqgMFAoAQDFTmqGVVEJqcq02uXlYValID/5eV9+VseJEydUrVo1TZw4UYMGDTI4KVAyUCgBAG6na9eu2rlzp44cOSKLxWJ0HMDtcZc3AMDtjBw5Un/88cc1n88NwPWYUAIA3FKzZs2UlpamrVu3ymS6+ULpAPKPCSUAwC2Fh4dr+/bt+vXXX42OArg9JpQAALdkt9tVp04dVatWTcuXLzc6DuDWmFACANyS2WzWiBEjtGLFCh0+fNjoOIBbY0IJAHBbGRkZqly5stq0aaNp06YZHQdwW0woAQBuy9vbW4MHD9bXX3+tM2fOGB0HcFsUSgCAWxswYIA8PT01adIko6MAbotCCQBwawEBAerXr58mT56sixcvGh0HcEsUSgCA2xs6dKguXLigmTNnGh0FcEvclAMAKBH+/ve/67ffftOxY8fk4eFhdBzArTChBACUCOHh4YqKitLSpUuNjgK4HSaUAIASo0WLFrpw4YK2b9/O4xgBF2JCCQAoMUaOHKkdO3bol19+MToK4FaYUAIASgyHw6H77rtP99xzj77//nuj4wBugwklAKDEMJlMCg8P1w8//KCDBw8aHQdwG0woAQAlSmZmpipXrqwnnnhCM2bMMDoO4BaYUAIAShQvLy8NGTJE33zzjU6fPm10HMAtUCgBACXOgAED5OPjo88//9zoKIBboFACAEqcMmXKqH///vryyy+VkpJidByg2KNQAgBKpCFDhujixYuaPn260VGAYo+bcgAAJVbPnj21ceNGRUZGytPT0+g4QLHFhBIAUGKFh4crOjpaixcvNjoKUKwxoQQAlGitW7dWfHy8du7cyeMYgXxiQgkAKNFGjhypXbt2acOGDUZHAYotJpQAgBLN4XCobt26uuOOO7Ry5Uqj4wDFEhNKAECJ9tfjGH/66Sft27fP6DhAscSEEgBQ4mVlZalKlSpq0aKFZs2aZXQcoNhhQgkAKPE8PT01dOhQzZ07V6dOnTI6DlDsUCgBAJDUr18/lSpVSp999pnRUYBih0IJAIAkf39/DRgwQFOmTNGFCxeMjgMUKxRKAAD+a8iQIUpPT9e0adOMjgIUK9yUAwDAFV544QWtW7dOx48f53GMwC1iQgkAwBVGjBihP//8UwsWLDA6ClBsMKEEAOAqTz31lGJjY7V7924exwjcAiaUAABcJTw8XHv37tXatWuNjgIUC0woAQC4isPh0EMPPaRy5cpp9erVRscBijwmlAAAXMVkMmnkyJFas2aNdu/ebXQcoMhjQgkAwDVkZWWpWrVqevTRRzV79myj4wBFGhNKAACuwdPTU8OGDdP8+fMVExNjdBygSKNQAgBwHX379tVtt93G4xiBm6BQAgBwHaVLl9bAgQM1depUJScnGx0HKLIolAAA3MDgwYN16dIlTZkyxegoQJHFTTkAANxEnz59tGrVKp04cUJeXl5GxwGKHCaUAADcRHh4uGJjYzVv3jyjowBFEhNKAABuwTPPPKOTJ09q7969PI4RuAoTSgAAbkF4eLj279+vVatWGR0FKHKYUAIAcAscDocefvhh+fv7a926dUbHAYoUJpQAANwCk8mk8PBwrV+/Xr///rvRcYAihQklAAC3yGq1qnr16mrUqJHmzp1rdBygyGBCCQDALfLw8NDw4cO1cOFCnTx50ug4QJFBoQQAIA969+4tf39/TZgwwegoQJFBoQQAIA9uu+02DRo0SNOmTVNSUpLRcYAigUIJAEAevfrqq8rKyuJxjMB/cVMOAAD50K9fP33//feKioqSt7e30XEAQzGhBAAgH0aMGKEzZ85ozpw5RkcBDMeEEgCAfGrbtq0iIyO1f/9+mc3MaFBy8bsfAIB8GjlypA4dOqSVK1caHQUwFBNKAADyyeFwqFGjRvLx8dHPP/9sdBzAMEwoAQDIp78ex/jLL78oIiLC6DiAYZhQAgDgBJvNpho1aqhevXpasGCB0XEAQzChBADACRaLRcOHD9fixYt14sQJo+MAhqBQAgDgpN69eysgIECffvqp0VEAQ1AoAQBwkq+vr15++WVNnz5dCQkJRscBCh2FEgAAF3j55Zdlt9v15ZdfGh0FKHTclAMAgIsMHDhQ3377raKiouTj42N0HKDQMKEEAMBFhg8frri4OM2ePdvoKEChYkIJAIALtW/fXocOHdLBgwd5HCNKDH6nAwDgQiNHjtSRI0f0/fffGx0FKDRMKAEAcLHGjRvLw8NDGzduNDoKUCiYUAIA4GIjR47Upk2btG3bNqOjAIWCCSUAAC5ms9lUq1Yt3XfffVq8eLHRcYACx4QSAAAX++txjEuXLtUff/yR/XpqhlUHYpO1KzpJB2KTlZphNTAl4DpMKAEAKADp6ekKCQnRU936qEqr57XhSJyiE9N05Tddk6SQQF81Cw1W9wYhql6htFFxAadQKAEAKAAxiWnqOv5bxTrKyGKSbDf4bmsxm2SzO9S0WjmNbV9HFQN9Cy8o4AIUSgAAXGx+RLTeWXFAVpv9hkXyahazSR5mk95rG6Zu9UMKLiDgYhRKAABcaNKGYxq/+qjTxwlvXUOvNKvugkRAweOmHAAAXGR+RLRLyqQkjV99VAsiol1yLKCgUSgBAHCBmMQ0vbPiQL7fn7Dyc50c94ziFr2X/drbKw4oJjHNFfGAAkWhBADABd5Ytk9We/6uIss4fUwX962TycMrx+tWu0NvLNvninhAgaJQAgDgpGNnU7QpMl62fBRKh8OhpLVT5Hdvc5l9y+bYZrM7tCkyXpFxKS5KChQMCiUAAE6asy1aFrNJkmTPytCpqQN1aupA2bMysvexpafoz4k9dWZ2uBx2W/brqfvXK/PcSQU8+vw1j20xm/TNVq6lRNFGoQQAwEkbjsRlTyfNnt4q98wwWZNidX7j19n7JK7+QvaMNAW1GSaT2SJJsmek6fzPs1SmURdZbgu45rFtdoc2HI0r+C8CcIKH0QEAACjOLmZYFX3VjTPed4bKv2FHXdi6RL41GsmWel5phzYqoEU/eQbelb1f8ub5Mnl4yb9+uxt+RnRCmlIzrPLz5ts2iiZ+ZwIA4ISTCam61pWTZZv8XemREUr4/lPZsy7Ju+K9Kl2vbfb2rMRTurBjhcq1HSmTh+cNP8MhKSohVWF3lnFteMBFOOUNAIATMq32a75usngq6OkhsiaflSMzXUFthspkMmVvT1w7Vd531ZRfzUec+hygKKBQAgDgBC+P638rvXTid0mSw5opa1Js9uvpUXt06fhO+ddrK+v5s9k/5LDJYc2Q9fxZ2TNynka/0ecARuOUNwAATqgU5CeTlOu0d2bcCZ3fPE9+dVoqK+6EElZO1J19Jsns4yfbhXOSpHPLxuY6ni0lQae+7KuAFv3kX/9ZSZLpv58DFFUUSgAAnODn7aGQQF+dvOLGHIfNqoQfJshyW5ACW/aXNfmsTn81TInrpqlcm6Hyuec+le8wOtexEn6aJA//YJVp3EWe5Stlvx4S5MsNOSjS+N0JAICTmoUGa/a2k9lLByVvWaDMs8dV4bn3Zfb2lVdwZZV95Dmd3zhbfjUfUamq9eVRJjjXcRLXTpPFr6x8azTKfs1iNqlZjdz7AkUJF2QAAOCk7g1CsstkxplIJf+2UKUfekY+99yXvY9/w07yuqO6ElZOlP3SxVs+ts3uUI+GIS7PDLiSyeFw5O/BowAAIFvP6du05XhCvh6/eD0Ws0mNqwRpdt8GLjsmUBCYUAIA4AJj29eRh9l08x3zwMNs0tj2dVx6TKAgUCgBAHCBioG+eq9tmEuPOaZtmCoG+rr0mEBBoFACAOAi3eqHKLx1DZcca2TrUHWtz7WTKB64hhIAABebHxGtd1YckNXuyNM1lRazSR5mk8a0DaNMolihUAIAUABiEtP0xrJ92hQZL4vZdMNi+df2ptXKaWz7OpzmRrFDoQQAoAAdO5uiOduiteFonKIT0nI8Uceky4uWN6sRrB4NQ1QtuLRRMQGnUCgBACgkqRlWRSWkKtNql5eHWZWC/HgCDtwChRIAAABO4S5vAAAAOIVCCQAAAKdQKAEAAOAUCiUAAACcQqEEAACAUyiUAAAAcAqFEgAAAE6hUAIAAMApFEoAAAA4hUIJAAAAp1AoAQAA4BQKJQAAAJxCoQQAAIBTKJQAAABwCoUSAAAATqFQAgAAwCkUSgAAADiFQgkAAACnUCgBAADgFAolAAAAnEKhBAAAgFMolAAAAHAKhRIAAABOoVACAADAKRRKAAAAOIVCCQAAAKdQKAEAAOAUCiUAAACcQqEEAACAUyiUAAAAcAqFEgAAAE6hUAIAAMApFEoAAAA4hUIJAAAAp1AoAQAA4BQKJQAAAJxCoQQAAIBT/h/Lz73mkf/zogAAAABJRU5ErkJggg==" + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "import networkx as nx\n", + "nx.draw(graph, with_labels=True)" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2023-10-31T09:16:33.274958224Z", + "start_time": "2023-10-31T09:16:33.258879935Z" + } + }, + "id": "bc41f2e657f4f715" + }, + { + "cell_type": "markdown", + "source": [ + "Calculating all marginal distributions is done by using the sum product algorithm." + ], + "metadata": { + "collapsed": false + }, + "id": "49cb617c5fd7f30e" + }, + { + "cell_type": "code", + "execution_count": 10, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "โ•’โ•โ•โ•โ•โ•โ•โ•คโ•โ•โ•โ•โ•โ•โ•โ•โ•โ••\n", + "โ”‚ x1 โ”‚ P โ”‚\n", + "โ•žโ•โ•โ•โ•โ•โ•โ•ชโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•ก\n", + "โ”‚ 0 โ”‚ 0.44998 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 1 โ”‚ 0.55002 โ”‚\n", + "โ•˜โ•โ•โ•โ•โ•โ•โ•งโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•›\n", + "โ•’โ•โ•โ•โ•โ•โ•โ•คโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ••\n", + "โ”‚ x2 โ”‚ P โ”‚\n", + "โ•žโ•โ•โ•โ•โ•โ•โ•ชโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•ก\n", + "โ”‚ 0 โ”‚ 0.324403 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 1 โ”‚ 0.169548 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 2 โ”‚ 0.506049 โ”‚\n", + "โ•˜โ•โ•โ•โ•โ•โ•โ•งโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•›\n", + "โ•’โ•โ•โ•โ•โ•โ•โ•คโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ••\n", + "โ”‚ x3 โ”‚ P โ”‚\n", + "โ•žโ•โ•โ•โ•โ•โ•โ•ชโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•ก\n", + "โ”‚ 0 โ”‚ 0.333492 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 1 โ”‚ 0.0419636 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 2 โ”‚ 0.297586 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 3 โ”‚ 0.326958 โ”‚\n", + "โ•˜โ•โ•โ•โ•โ•โ•โ•งโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•›\n", + "โ•’โ•โ•โ•โ•โ•โ•โ•คโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ••\n", + "โ”‚ x4 โ”‚ P โ”‚\n", + "โ•žโ•โ•โ•โ•โ•โ•โ•ชโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•ก\n", + "โ”‚ 0 โ”‚ 0.211013 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 1 โ”‚ 0.165152 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 2 โ”‚ 0.198958 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 3 โ”‚ 0.14797 โ”‚\n", + "โ”œโ”€โ”€โ”€โ”€โ”€โ”€โ”ผโ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”€โ”ค\n", + "โ”‚ 4 โ”‚ 0.276907 โ”‚\n", + "โ•˜โ•โ•โ•โ•โ•โ•โ•งโ•โ•โ•โ•โ•โ•โ•โ•โ•โ•โ•›\n" + ] + } + ], + "source": [ + "graph.sum_product()\n", + "for variable in graph.variables:\n", + " print(graph.belief(variable).to_tabulate())" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2023-10-31T09:16:33.319109726Z", + "start_time": "2023-10-31T09:16:33.277484723Z" + } + }, + "id": "47a7412bada04433" + }, + { + "cell_type": "markdown", + "source": [ + "The joint most probable state of the graph is calculated by using the max product algorithm.\n", + "The result is a random event that describes the most probable state." + ], + "metadata": { + "collapsed": false + }, + "id": "ce9779d3638d933a" + }, + { + "cell_type": "code", + "execution_count": 11, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{Symbolic(name='x2'): (2,), Symbolic(name='x1'): (1,), Symbolic(name='x4'): (4,), Symbolic(name='x3'): (0,)}\n" + ] + } + ], + "source": [ + "graph.reset()\n", + "print(graph.max_product())" + ], + "metadata": { + "collapsed": false, + "ExecuteTime": { + "end_time": "2023-10-31T09:16:33.319755750Z", + "start_time": "2023-10-31T09:16:33.318917680Z" + } + }, + "id": "72f3105bbcd8d99a" + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 2 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython2", + "version": "2.7.6" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/requirements.txt b/requirements.txt index f721b39..4fea548 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,4 +1,4 @@ -networkx==3.2 -numpy>=1.26.1 +networkx==3.1 +numpy==1.24.4 random_events>=1.1.1 tabulate>=0.9.0 diff --git a/src/fglib2/graphs.py b/src/fglib2/graphs.py index 3776333..729bcef 100644 --- a/src/fglib2/graphs.py +++ b/src/fglib2/graphs.py @@ -439,7 +439,7 @@ def to_latex_equation(self) -> str: return r"P({}) = {}".format(", ".join(tuple(variable.name for variable in self.variables)), r" \cdot ".join([str(factor) for factor in self.factor_nodes])) - def brute_force_joint_distribution(self) -> Tuple[np.ndarray[np.ndarray[int]], np.ndarray[float]]: + def brute_force_joint_distribution(self) -> Multinomial: """ Compute the joint distribution of the factor graph by brute force. @@ -451,15 +451,16 @@ def brute_force_joint_distribution(self) -> Tuple[np.ndarray[np.ndarray[int]], n """ worlds = list(itertools.product(*[variable.domain for variable in self.variables])) worlds = np.array(worlds) - potentials = np.ones(len(worlds)) + potentials = np.zeros(tuple(len(variable.domain) for variable in self.variables)) for idx, world in enumerate(worlds): - + potential = 1. for factor in self.factor_nodes: indices = [self.variables.index(variable) for variable in factor.variables] - potentials[idx] *= factor.distribution.likelihood(world[indices]) + potential *= factor.distribution.likelihood(world[indices]) + potentials[tuple(world)] = potential - return worlds, potentials + return Multinomial(self.variables, potentials) def reset(self): """ diff --git a/test/test_graphs.py b/test/test_graphs.py index 5ad2d70..60fa6cd 100755 --- a/test/test_graphs.py +++ b/test/test_graphs.py @@ -123,11 +123,8 @@ def test_graph(self): self.assertEqual(len(self.fglib_graph.edges), len(self.graph.edges)) def test_brute_force(self): - worlds, potentials = self.graph.brute_force_joint_distribution() - for index, variable in enumerate(self.graph.variables): - for value in variable.domain: - indices = np.where(worlds[:, index] == value)[0] - print("P({} = {}) = {}".format(variable.name, value, np.sum(potentials[indices]) / np.sum(potentials))) + distribution = self.graph.brute_force_joint_distribution() + print(distribution.to_tabulate()) def test_calculation_by_hand(self): x1_to_fa = self.graph.node_of(self.x1).unity()