templates/bsda_colab_template.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "bsda_colab_template.ipynb",
      "provenance": [],
      "collapsed_sections": [
        "Z7_dLXw_yM9e",
        "GMlmnXt6zQjC",
        "xHC-Y7Rryx_H"
      ],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "ir",
      "display_name": "R"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/MansMeg/BSDA/blob/main/templates/bsda_colab_template.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# BSDA: Assignment COLAB Template\n",
        "\n",
        "This template contain information on using Google COLAB to run code in the assignments of the course. Just as R-markdown, COLAB is built upon markdown. More information on how you can add equations and format your report you can find [here](https://colab.research.google.com/notebooks/markdown_guide.ipynb).\n",
        "\n",
        "When using this template for assignments report, just keep the \"Setting the stage\" part below. That will install everything you need. Then start by filling out the general information. After the general information below are examples of using `bsda`, `markmyassignment` and `rstan`. Remove these parts and add your own work instead. Good luck!"
      ],
      "metadata": {
        "id": "uYMilov8ubkU"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Setting the stage\n",
        "We start by installing nessecary packages for the course. This include Stan/Rstan, loo, the bsda course R package, and markmyassignment."
      ],
      "metadata": {
        "id": "Z7_dLXw_yM9e"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### Installing Stan/Rstan"
      ],
      "metadata": {
        "id": "sDU8qkbfNKLB"
      }
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "rCEtur7fajGn"
      },
      "source": [
        "#cat(system('python3 -c \"from google.colab import drive\\ndrive.mount()\"', intern=TRUE), sep='\\n', wait=TRUE)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3o11KxDzyC0g"
      },
      "source": [
        "# takes rougly 1m \n",
        "#install.packages(\"rgl\", repos = \"http://cran.rstudio.com/\")\n",
        "#install.packages(\"ConsRank\", repos = \"http://cran.rstudio.com/\")\n",
        "#library(\"ConsRank\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "SPYGePw53fWq"
      },
      "source": [
        "# system(\"add-apt-repository -y ppa:marutter/rrutter\")\n",
        "# system(\"add-apt-repository -y ppa:marutter/c2d4u\")\n",
        "system(\"apt-get update\")\n",
        "system(\"apt install -y r-cran-rstan\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "1IzjHJFxv01f"
      },
      "source": [
        "#install.packages('ggmap')\n",
        "#install.packages('codetools')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "sdyV67dk2vGx"
      },
      "source": [
        "#library(\"ggmap\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "qIbhtDzdwYdB",
        "outputId": "fa4dc4bc-93e2-4fb4-ad09-553d8f0ee2e8",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "source": [
        "options(repr.plot.width  = 3,\n",
        "        repr.plot.height = 3)\n",
        "Sys.setenv(USE_CXX14 = 1)\n",
        "library(\"rstan\") # observe startup messages\n",
        "options(mc.cores = parallel::detectCores())\n",
        "rstan_options(auto_write = TRUE)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "Loading required package: StanHeaders\n",
            "\n",
            "Loading required package: ggplot2\n",
            "\n",
            "rstan (Version 2.21.5, GitRev: 2e1f913d3ca3)\n",
            "\n",
            "For execution on a local, multicore CPU with excess RAM we recommend calling\n",
            "options(mc.cores = parallel::detectCores()).\n",
            "To avoid recompilation of unchanged Stan programs, we recommend calling\n",
            "rstan_options(auto_write = TRUE)\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "### Installing `bsda`"
      ],
      "metadata": {
        "id": "mcJ1FPcvNOS4"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "install.packages(\"remotes\")\n",
        "remotes::install_github(\"MansMeg/BSDA\", subdir = \"rpackage\")"
      ],
      "metadata": {
        "id": "UvPhBYf9NSo6",
        "outputId": "27e00a3a-58e1-487a-b58e-8f695ddcb9c4",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "Installing package into ‘/usr/local/lib/R/site-library’\n",
            "(as ‘lib’ is unspecified)\n",
            "\n",
            "Downloading GitHub repo MansMeg/BSDA@HEAD\n",
            "\n"
          ]
        },
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "\n",
            "* checking for file ‘/tmp/RtmpZR0NPV/remotes3d54992b8c/MansMeg-BSDA-53218b3/rpackage/DESCRIPTION’ ... OK\n",
            "* preparing ‘bsda’:\n",
            "* checking DESCRIPTION meta-information ... OK\n",
            "* checking for LF line-endings in source and make files and shell scripts\n",
            "* checking for empty or unneeded directories\n",
            "* building ‘bsda_0.3.1.tar.gz’\n",
            "\n"
          ]
        },
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "Installing package into ‘/usr/local/lib/R/site-library’\n",
            "(as ‘lib’ is unspecified)\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "library(bsda)\n",
        "data(algae)\n",
        "head(algae)"
      ],
      "metadata": {
        "id": "FPCSB9tlNui2",
        "outputId": "b72f0916-82d7-420d-91cd-c9c8287366ed",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/html": [
              "<style>\n",
              ".list-inline {list-style: none; margin:0; padding: 0}\n",
              ".list-inline>li {display: inline-block}\n",
              ".list-inline>li:not(:last-child)::after {content: \"\\00b7\"; padding: 0 .5ex}\n",
              "</style>\n",
              "<ol class=list-inline><li>0</li><li>1</li><li>1</li><li>0</li><li>0</li><li>0</li></ol>\n"
            ],
            "text/markdown": "1. 0\n2. 1\n3. 1\n4. 0\n5. 0\n6. 0\n\n\n",
            "text/latex": "\\begin{enumerate*}\n\\item 0\n\\item 1\n\\item 1\n\\item 0\n\\item 0\n\\item 0\n\\end{enumerate*}\n",
            "text/plain": [
              "[1] 0 1 1 0 0 0"
            ]
          },
          "metadata": {}
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "### Installing `markmyassignment` and `loo`"
      ],
      "metadata": {
        "id": "2dybTe54N9Hi"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "install.packages(c(\"loo\", \"markmyassignment\"))"
      ],
      "metadata": {
        "id": "AHqaB1aCOD2d",
        "outputId": "a2c83752-5847-414e-82e3-4d34b4a55a41",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "Installing packages into ‘/usr/local/lib/R/site-library’\n",
            "(as ‘lib’ is unspecified)\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## General information\n",
        "When you write your assignment report, don't forget to include the general infromations. See the R-markdown template for details [here](https://github.com/MansMeg/BSDA/blob/main/templates/assignment_template.rmd)."
      ],
      "metadata": {
        "id": "d5lVU35CUGZU"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Running R code in COLAB\n",
        "We can run R code in COLAB just as we would run R code in R-Studio. Just add a code block with your R code and press the run button."
      ],
      "metadata": {
        "id": "GMlmnXt6zQjC"
      }
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "xLlQFtlPxD5j",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "46df5bff-7aff-43f2-90bd-97a366653934"
      },
      "source": [
        "x <- seq(0, 2*pi, length.out=50)\n",
        "plot(x, sin(x))"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "plot without title"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAFoCAMAAABNO5HnAAADAFBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUW\nFhYXFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUmJiYnJyco\nKCgpKSkqKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6\nOjo7Ozs8PDw9PT0+Pj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tM\nTExNTU1OTk5PT09QUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1e\nXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29w\ncHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGC\ngoKDg4OEhISFhYWGhoaHh4eIiIiJiYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OU\nlJSVlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWm\npqanp6eoqKipqamqqqqrq6usrKytra2urq6vr6+wsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4\nuLi5ubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnK\nysrLy8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc\n3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u\n7u7v7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7////i\nsF19AAAACXBIWXMAABJ0AAASdAHeZh94AAAbZklEQVR4nO2deWBM1/7Av3fWTDJJZJUQERE8\nS2JpbI2IrU0tpfTV2qIJRXi0aKhayvMeDVVLeKoPXWxtKaq0r1o0Xqs/RFFbLUURVNDYItuc\n38zcmWSSzL1z79yb78zknc8fc4659545PjO595xzz/leIBQUwNUV+F+BikaCikaCikaCikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikZCimjDhd2ff/7d77LVpVrjvOg7\nk0LBTOScRzJWqLritOicetBg+KyMjOmDakHzO3JWqXritOhU9aeWXPFyZoJMtanGOC06LKUs\nP6COHFWp3jgtWv2PsvxbGjmqUr1xWnTd/mX5PlFyVKV647ToCcyCx2zuwUyYwrfn0cPVjqPi\nfTkt+m4r8O06fNzYYZ28IfE+z46HoBpyCE80KVjUQmn6THW7VcV8+/0ABSJLvjTnrwk1vNok\nagFU0PVXp2tYZRTAD6KPkdQFzz+bnX3OkUbRotfpmg/VhivGxjUdodqWEaoI7PIv3m8SH3TR\nghAn+vtB0dBkdUpCyfeawDtkUqB2ZLOW6QHd8quses7gRqJzLpTymRjRs5UDmseP99e+R0io\nbzFZCFvIfuXlxb7xq65XTUWdwlWiz3ftWvEdxvbKIXwo5EvVlyR4I7kILxCDCk6RtvATyYOA\nkDh1Hd0qGWoqE64S/TNUKuX3sl/0XOBrk9hgWB2vYGKmaL4hxM8rn/jAAaKFc2QfvFyQAydX\nqrbJUFV5cJXo/F9+4dm6UqBow1D9dN+/L4tRLyDkeThCWsOVEkWYgbRRF5A7cCx/WL17MtRV\nFtzoHG2DUNHrvH8m3jtJXkhAAcmGoUVPMbs3MYPJecVzhOxRD1EBMD0vVm1VheIy0bnneDYK\nFF0cP/ghiZ1PyCbmqUuHGH+9oh7DNA4c6Qd7SGGCX+yuS7CmS82LclRXMi4TPYWvFEGiH77m\nA6DoPjEilzyEZqDzVgX0Th3Vs6HGZ5z/2j1JPpF3yS9wrbBjPzmqKxnPFV3QIWqT38aspwOa\nNPnqD/igs25mFtsmzBvtZfwC+kRkEjIxlpDdarc4T3uu6HdDrpF2b5CS7p1S1MaWYccTZZsK\nd3il5TH7S5aqvtzQLRhaZbpBLxFV9BM2hEkVHf8WIe/5nyWHmJzc2MRb5TfuCQtgWtfyWTvc\ne8JyGBPoBr1EVNEKhbYUpTTRJVm6F7flFT8bvPj/mLdaRl6puP3hp7Vbf3Drfd8j5N/BRZdr\nT3O2yrKBKnqKb1lTQ9qp41RztaqBf42Pi96uB+Cb+oedXdbr9pMW08m5sJmErA4qcqrCMoIq\nurBlfKE1L0n0zfA+N5JHFixUbSVkg5Zj3/GalxXTJvj1Ml4iL8EFZ+orJ7gXw1O6ydasJNGT\n4grJDvW3ZFq04XrMaK69vuwNoclrSsjBtLbwwrfiaysryK2OvNvW3L55PLs5Ev2XxcaXqarU\nd2B0aHue1lvUMuPLXGX3F1V9VGkGcVWVGc/sgvubB4u+6h0FcYv4BlRnROWSbert+W2GkB99\nl8lZQ9F4oujdfVX6JLPgW3CEt5z7LRt91nrEzraR1wjJiHTpT9oDRc9SDevQ7M3QdnmELAt2\n0JrIG+sNoBmUY8yegUptQEw8T/R/VLvIGd3UnEYjyDf6xQ6Lugcb2ZbOdTgjW/2cwPNE9xpm\nfPk6sG6SIp6ZIuBsEPKh6fWbgfWYPptlqp4zeJ7oWutMr7czU2DkMSFljWn5mJCJqsFNY1N1\nQ1w36OF5ooM/Y9MS1R5BZV2P6Hp8nW5Df/9T5FjgAjlq5xSeJjrv/Zpx80+acsdA4LqB37qA\n0gueME3JWlzbZS0PDxO9NzS8hbaZYprxB92no+DizjNzzpozv7qu5eFZoi/oxxcWJtV9wzvj\nh541TtjdxR4P4CCbuQpn5aieM3iW6FEJxj/9RxN9AJinTosor+Yac3J/qfqGDLVzCs8S3WC5\nOSk6CN+LKm9804eE3ByoAIZ59rLk2jmFZ4kOsraEtf8RVd6t+u2ycmIax4eczkqq7ZrTtGeJ\nbrqITa+DoCZ0GTn9FAoGnjpv/A+3GyKxcs7hSaKvv93Yd1KWKTe7TonYIu/WSGfvK36heyyt\ncs7hQaJ3+DV6ySdUMaK4eLlqk+gi8yCbzVyC3yTWzik8R/QZr5kGcqKZRhUVqH9ffJGPmf+y\nmdOQI7F2TuE5olM7m15LvhukWX+70kYBNJ/Bpu+4pnfoPqJ//0t0KcFQ+f5UgxVs+oARX2MT\nq/U/EVL4djQoGs52wTQP9xFd8MF7pQyx84sO/cSS0e1yqnzDaG3qiiY6dZf979Zu/0BCRZ3D\nfUTbYu/U0dJyN/ea2KZdKdv7Bim7bjSeOG5EpTtfNyfxGNEPX6nFzs2YWs/5c2zUu2y6Jhh9\nYNpDRN8fowEGYg+RR3NV250u9oF1UeVZuCqhdk7hGaILE2K+uHchHpj6mpBPOY4RwD04zGbO\n4Q+XeobozGBz2/doQsS3D6WUG5nJph8HoM/F8wzRHd5g0/MS72S/FWFaemjICnvmZ+zGtGeI\nrrXekvF2rmln5VFCxIojmyNBXQtiD0sqSTSeITra0ucuVAq7IctJ/oxIAK9RD8i1If5i7hxI\nR7roi7s3rty4+6JM9WGpKHrg82y6Sy095NULCeamnaF7b8lFiUGi6At/i7GsKY4ZL+OgWAXR\nx2cx00wXwWsxIyQXXeKzhc18rUHtiEsS/UeKCkIHTF/60dLpA0JAlWJv4r1TlBP9e2eoEwaa\noWtfDeSNpiKMXLAs2f0NLkkuTARSRO8NVQ7Jtl69DdlDlKH7ZKqVreh7DRLPErK/FRPaY6UM\njbLHTBabyQbUyHtSRKs7ni/3/rmOahlqZMJW9NwodgRoeqTouyp2afMam05vIktxQpEiekbF\nAYPiGZLrw2IruvUcNr1hvUUika1q00nakKkefQSzLS211bHfmilZJEd1LNiKLm1D+34hT+Hz\nlZ2mpvqCXyQ0k+erE4RU0YrX2BAmZxPkbF3bim68lE0fKcTN5eDm6OREdSNjzydnsD/elGmp\nortDowPGn/NinXaOfJUqJ3q0ZY7dem9Jwxzl6N/FfL43JPeVrUhHSO6wbK6jeP2XROgscFKb\nsPjRNqJvf6wZY2psHA6eLraenBTpLEOtu7Riw745jfSe4YPJKghaK+xIofGjS0X/OUyp1oOm\nV3qycrh8Y/U34SSbQRwulS768Swt1BQ2z0Jw/Gir6IK2jb8rJEeeVLSdvFd0Nbl5ZL2/exDy\nZCyWF8midzdg0rLbQ08hcwcFx4+2il4WctOcjmskb0ssfhKbTo2TtVg+pIoeBFF7jBfDBV4+\n7zo+UHD8aKvojpZIsVfBiRCqPGxRf25KtmoyPOYczYxmlZxpL6B5Jzh+tFV05IeWN/xljqs2\nT5mUPj7MeLHQvCTbAA0/UkWXrmUvWej4QMHxo62iG1tuPRVpdouonxCOpXf1Dl2Sm7uzVcwt\nx3vLgBTRle9ROLhrITh+tFX0y8+w/94pwzh0RcY2NX/Gg2ZjZC/aHlJEa5ZW2LDUQTh5wfGj\nraKzVeaf9KV6o0RX0hElAevYzIYaVTvH4za7cFeK6HHQ3vYW3q72MNbBkULjR7Oi17ZQAxOW\nvmSEbzf5OoVWStvSp6EKY8VeT60J6hZricRz9IYgaDxh24kb+TdObJvQGII2CDhYUPxos+jR\n3jP3HssI1TTq/5E8A6TluGsNjPBzFY5Lnw+PX3ds70zvMVIvhvfn1SyNjxs2X/rtDysm0Ts0\nB0zZwm5dZCu2HNEZbLqwXtWUb6JLN/N544D6S8mtjpKD81J6d+qdMu+QnL86k+jnhrP5Y1UU\nD2lRoPnccSpQQHPJSS7AcTYz/Dk3ih99J+2VUhKNohtYgz/77pDh8ypT3F//6qaNyUqvsG4f\nVdEdgC98LZlVDdwofnQl0db1E1Ukmhg+7BbmpUzeuOlVff+qaXjIKdrwaa8WTVlEleE4fnTZ\nqaPqFvcsDDplSk4Eynl7qAw5Tx0LALz9WeSpnBmMi6GJKMv5OSO6asqX8WIYkSzuUiVi4N/c\nvFsXH151sQFvW8eqsuFP+Usv+ah/I01oxnE5mneEqH8Sc6S4gX9Th6VmalV2J8Ay/+4k3JS9\n8IfdfEcseT2MkaPDYiTigIgDRQ/8Fzq10E0wxTU2spl1gfJ3iUZFm2dCZarYW+1SRb+eJuJA\n0QP/Vc2YWPPknPtNHY0diOe2eiebSX7ZnEgVfT958NenzplxfKDogf+q5lZMq525uUtCfbqm\nO7vUi4tvNJb5a8samxOpom2eUeP4QNED/1XOHy9pTKMH49OTlHyxUp1gWw1L5sO65kTyraxh\nqVYcHyh64B+BggzNVlNqubslG0eZa2xmShL7OZg9Q9ED/xjETmXTifGyFmtoOM6c3ghhg89K\nEX3d2HC4XobjA0UP/CPwpzWw1X8ZOR9VvndyW8WTR0jhd43bsmPCUkRDsshztMiBfxSugOUq\nLmdbuni4Mjm9lwb0auUwS1dIiugBxuvHgDKEHSx84B+Hx9qv2Mw2b/lGlqaHmO6eFo3RrCvt\nCHjGqqyqpO8z5kHSki4CfykCeOBtmWrcsexep3TRpt/B45/kndSNKfq0/5AcQnY1VCdOlmuW\nzj6lZR3SkrJVBVJFF6f9lZCL0QAd5FSDKZpkN4NIP/BNndpJOV+eEr/ws2TW1yp9T6roeTCR\nkB7MmDSFnA1+VNHEcGS0OtP0F7lFvVWWArMZy3V1dpvS96SKbtaPkKuMsbOS0kJi7WzBFU1I\nE8vc61fb8O8nkJJItrz7UXNL35MqWr+SkNXwLSHLa3DuLx5k0betUw++Z2SIiVe0skco02o/\nIWcTG5SFhpIq2tcoepCPsbmW6SO9iqUgiy4NhXccciUXdj8x8NW1QzUQVge62NzgkHzqGExu\n6E0rQUY2klzFMpBF52ssoU636KWPS49oYBrkeDiNmXXc9m2pov8J7WvBPkI+1LwusYK2YJ+j\ne3c3t06Ln+zvaE+H3FFbpsn1G1jufami84fr/E1zHcNj5ZxYhS36lN+Qa+T+KC+AOjMlLsXf\no7KMQ6+qX+59uXqGB2SNnYMtmhyOhdpqiNx8ZEVEgrSxpZ3eloxNG9oE7YKzGH5ODssynT9y\nar8lqaAz1hlsUxPLvU9FWygK+JjNLIuUUszDbyMSzD35nKDMchuoaAu/gyVUwyEpn/1piKY+\nA/EX7n1Rv0NhuS1UtIWr1nHpg+B84NLtqrmPyMFYYECTVuF/QEVbKA5ey2YWOT9f2lCPDRt3\nodYrlZYoUNFW0qNuGFVt6KoM6utsUM6j1hjg/2hVaZv7iHYU1riqedCu9rv7O2t0TVakaEY7\nN7q+y9q0+yS00jb3Ee0oUHeVkz+7oQKiM4zXsAM+q50q4QfGcsZY3rDSNvcRbYsrTh1GDLXe\nYTPTmztzeO567aDvzIMlnSuHjaOibbhmnV26n3Fijfj7+sAolbrZSWKY4VU5sg0VbcNvYInR\ncNiJU9dG1fLi4hHKUJ+hjfztLAahom3I17GGbqUHiB4iK6ljjnWUNcm38dv25hJR0bYMaW88\nZZzrBoxC0U/k01qOWhfgvmN/MQ8VbcuVWkn7T4XEN47JyWpbX9zNlq+1lszmILvbqehyXH6W\nYUAx8CYhD5uOF3XkPjjIjhRnVm7amaCiK3BDvZT97DU1RRx1uhsD4DPxkbGJmDDa7h5UdAVK\nH2ghZnDplxq9fsjwfqNuUmHBeL39ZZFUdAWuWB5Le3YOI2CxiIXEPsZ+yjRFM22L8NC99neh\noitgqL3E+Hr0CfBSQleBa3Yvsw87Ojk/Lux9rvhuVHRFMgKPkVP+Azbo1h3vGiFskeN3Ksss\nhc+COfehoitSPFiXGtt0sGoiIY9bCoq3dGwkpGSaJ0J/XJtzJyq6Mpv7MPUGfmPKfVB5uLMS\nhilMvCKpbuDXxvywZzl3wxYtOgisK7CGWSrcKKAei/XfkNRG16bofiW7VNwPd8YVLToIrGu4\nAr8aX3MGagC8xzqIXloUnElIXtuabzbrMFQ1i3s/VNGi14K7CEPkAkKuRbbdOaL1p41a8Vcl\n2zwvsmBRkl7Vly8KIqpot1sLzsVSvx/JkDb529XbSW5dvgcXFCyKg6hnzcuNtvEHLEEV7XZr\nwbkwjFE/p3rxGaVpSfVSnrvi99qHjoZFKSrTze/Ff+EtElW0+60F5+TbF6Bt2kFSsiY5FHp/\nybXXqAY3DNHTyG71DlIYO4m3QFTR7rgWnIsLprstBb38JryheFnDMZB3T7uRkK2qhQUjk2/0\nDudfD0rXgnNQFLSGkBnh58j0lmS/zl6MypspvgD13i76uIZ/A5W6+Sn+8lBFu+NacE6mRVwu\nCvo3ydb/+9YHrWt/Umnq0ZXIVrOYg4uDny3O2/aiLsvRUgHcdrQbrgXnJL9bYBqsmKAbvtan\nVmtFQHiFJyneSoy7ncMcJmf9VxEyy3FkBPQuuLutBeemOLMVBHfbsEO1tCSLyUvT7bCZHnoi\nEYBRpyT1KCFT25NrIY7j7tOxDh7uqb8lpJmxNZFZp7cCwGsM20ssyJqp67yWufVVkyYByVkb\n/TZFJTqeBkJF89EvqegqnCB3I3067V0TMM5H07Dn6NTeASpvHTSDRyQ3YmJ34xfg85qA8Nbu\nEwTWFncRfbFmlzVwaVesX4dCsgf8RgY2Zpgohd6fOXzpKeYTQubFkoeDWwsKP+E+QWBtcRfR\n5LeeDIBqsNp4JeyjuEAGM5t2M08XvQQ/k4KAkDyy04cc8V4nqCRXiXYcBNZNuFdvWP5RuEMK\nVG2Ioaai5HdoTY7A84RkqGOW/d3vTf1QYVN83egc/fuFUua6j2hjz2/lCbj+eDh8T86BlhwA\nPcn38jMFtpsSwyjiVwucSu0y0bkV7zGfZ2ziM4GcgaQkskpXRx0XEgj3yE/QjpwElYG8ABfJ\nxmDypWqn4FJcJnpKpVJyyn7RnwHao5UEcH1VvO/iy8r/ks3wDin2DSXkPa3/+Cda9FeKeIaj\nG4m24Qe3Em3uJaa3jB2pCZxM7gRovi9JSF3dBKIHiXlkKBUtiOJ/dQlUhGZsVY1sGpemCNcO\nba4T1tgoBVX0EzaEeZRoE792BRWAtmMbrxoJf50j9untqKIVCm0pSr5SDkE15BCe6Cm+ZU0N\n3lMHOXrYTHr4x3ykq3k3Z0Am3+Y1MIv3cL+xvJvrD+Dd3K7P4fI4EevNadGFLeNLR8D4RVtY\nyR/WpnTyt30cPAUrH/ijuTt4vmtb/thtA2R4qJzzF8NTusnWLBXtGAmtjrzSWJ37hMTJo6KR\noKKRoKIlsSBB6J5UtCRGCS6AipYEFS0QKtoEFV3K/7zou1eE7rkmlnfzHl/ezRcY3lXGhcps\n3sPr8D+YJZH/GYgv8U1LFghe866A/ysxOFie5uABfQ42X+IPTpnDfxMoV4YH9+GJ/h+HikaC\nikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikaCikYCS/TdCXXV4ancUW8L\npyqe4D76zqRITVQfzrtVF0ZGa4L7/B/f578GqRxb1lomiP6d69BdHfX+nffyFS4EJNEFreD5\nf6So63HFBDrVypdH9O0o6DljiMrruP3NZ4I0L84aolb/yP35h5Scot+FQVNM7OHYvgbqT58c\nohE/I7o8SKIXwdvG10+AI+BIni7+nJZb9FhYZnzdAj3sb36KMa2L+By4n6tX1KI5p+hZ/JOd\nb+pbPiDknD6NbycBIIlu4WsO7hETan+B2e1JhYRH9KtdTTOEDbq69jdPNz+apljN/dSE+cxX\nnKInAG/c0gVgCoJHJD/BG0d0vpJdxDwcuG+i8ohmeazmnX52FZ7j2nReN+Yup+hhcKv4yi3O\nYpN1heSxgzBuQsARfRaGm9NZwB3lzKHoJeYTCAcP98b5cp4Cuob/yS36OXgzAKDheo7NdZsc\nSWCg/loHlXMIjuhsGGtOFwD3/ApHovdpOnDPGfAHeJHzj2UtbCbcojtB9LyP3vCDlfY3+9YN\nn7R5SSRwfRFCwRI9zpxmAPdT0R2I3qBtdZt769RXnlR04DB9M7AX4RH93WZTEO+T2kD7i8e0\n8KHxNUcfJijuATc4os/BMHM63fTUcQ54RRtmwjMOHqey1yfOfiykgfrLfKIt9IWDdt8PUpoD\neLwAHG1LoeCILlB1MqeDrE+hsQOfaEMK/M3hL2ow2A3utQtmXLly5SQMusJ7SRsF9hvSTyjN\na6LSnFhaWA6k5l1bb9PvoqQWT8hHPtET4J/cG6/GvWRO+9lvEE8qXRxoN2jc/RVseLYOHA2i\ncfCTKXkaHEUVdgCS6FVgeszuv2A29y48orcA7yzDCI3Jxa96vd2nU5/aYWITPL3jtL3NJbX1\npve3QUv7hR9muhi7AIcUcXw1EACS6OJE6DN7IBPLEbBon7ELrAwzvtifMlof/mbuJU+x34Pf\nqlQPfHO4D2Ta3crCfY7ezvikzujL+HFNR30VWsweqdPs5SlcCFiDSvcn11XXHsvVbphn/eu2\n30kr/eO/aP/wn54LUdbo9gXf5/NcDH/sXkNVayhn99CwsrmXfw/7V0oR0GFSJKhoJKhoJKho\nJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKho\nJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKhoJKqL6N3MIFPS\nXbHf1TXhoLqIJqNNESo2w2uurgcX1Ub0/agGjx/UaehGzyAvT7URTfYwsyYreOJXuZjqI5qk\nadWvu7oO3FQj0dkAv7i6DtxUH9El7WsGJUoO5lVlVB/RC2DTWljs6lpwUm1E/6rrQUhn77Ou\nrgcX1UV0SXufS0bb2gT7cZVcT3URvRDMT8qbA++4uiYcVBfRbg8VjQQVjQQVjQQVjQQVjQQV\njQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQVjQQV\njQQVjQQVjQQVjQQVjQQVjQQVjcT/Ax1yi047yrcEAAAAAElFTkSuQmCC"
          },
          "metadata": {
            "image/png": {
              "width": 180,
              "height": 180
            }
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Using the `bsda` R package\n",
        "When `bsda` is installed, it is just like using any other R package. "
      ],
      "metadata": {
        "id": "xHC-Y7Rryx_H"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "library(bsda)\n",
        "data(algae)\n",
        "print(algae)"
      ],
      "metadata": {
        "id": "qKMwKjjQO2Xn",
        "outputId": "715f1f7b-db0d-4e33-e61a-5d0b3ad07f76",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "  [1] 0 1 1 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0\n",
            " [38] 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1\n",
            " [75] 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 1 0 1 0 0 0 1 0 0 0\n",
            "[112] 0 0 1 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0\n",
            "[149] 0 0 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0\n",
            "[186] 1 1 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 1 1\n",
            "[223] 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0\n",
            "[260] 0 0 1 0 0 0 1 1 0 0 0 0 0 0 0\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Using `markmyassignment`\n",
        "The `markmyassignment` R package is a way to use use unit tests to check your code for automatic feedback. This will help testing that your code works as expected. For detailed information on how to use `markmyassignment` as a student, see [this introduction](https://cran.r-project.org/web/packages/markmyassignment/vignettes/markmyassignment.html)."
      ],
      "metadata": {
        "id": "utmepSbnOfL4"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "library(markmyassignment)\n",
        "assignment_path <- \n",
        " file.path(system.file(package = \"markmyassignment\"), \"extdata\", \"example_assignment01.yml\")\n",
        "set_assignment(assignment_path)"
      ],
      "metadata": {
        "id": "Eflx3mp5P69e",
        "outputId": "8eba450c-66ec-4ed7-f972-59e32c90eb05",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "Assignment set:\n",
            "Test assignment 01: An example of an assignment.\n",
            "The assignment contain the following (2) tasks:\n",
            "- task1\n",
            "- task2\n",
            "\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "In this example assignment, there are two tasks and also a mandatory requirement.\n",
        "\n",
        "**Mandatory requirement:**\n",
        "\n",
        "Store your name in the variable my_name.\n",
        "\n",
        "**task1:**\n",
        "\n",
        "Create a vector containing the values of π and e. The name of the vector should be task1.\n",
        "\n",
        "**task2:**\n",
        "\n",
        "Create a function that takes a numeric vector as the argument and returns the sum of the first and last element. Name the function task2.\n",
        "\n",
        "We start to solve this assignment by solving the first task."
      ],
      "metadata": {
        "id": "Vv29gtvKQTa4"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "task1 <- c(pi, exp(1))\n",
        "print(task1)"
      ],
      "metadata": {
        "id": "xsW0_S8WQLr_",
        "outputId": "286c3166-0670-4670-d380-4427a002d583",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "[1] 3.141593 2.718282\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "mark_my_assignment()"
      ],
      "metadata": {
        "id": "6a47bnKtQdwl",
        "outputId": "26f3018a-2716-44ff-9b70-c560e7df8eeb",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "\u001b[32m✔\u001b[39m | \u001b[33mF\u001b[39m \u001b[35mW\u001b[39m \u001b[34mS\u001b[39m \u001b[32m OK\u001b[39m | Context\n",
            "\n",
            "⠏ |         0 | mandatory-1                                                     \n",
            "⠏ |         0 | Mandatory tests                                                 \n",
            "\u001b[33m⠋\u001b[39m | \u001b[33m1\u001b[39m       0 | Mandatory tests                                                 \n",
            "\u001b[31m✖\u001b[39m | \u001b[33m1\u001b[39m       0 | Mandatory tests\u001b[90m [0.3s]\u001b[39m\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\u001b[1m\u001b[33mFailure\u001b[39m (test-mandatory-1.R:6:3): Mandatory tests\u001b[22m\n",
            "exists(\"my_name\") is not TRUE\n",
            "\n",
            "`actual`:   \u001b[32mFALSE\u001b[39m\n",
            "`expected`: \u001b[32mTRUE\u001b[39m \n",
            "Variable my_name is missing\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\n",
            "⠏ |         0 | task-1-subtask-1-tests                                          \n",
            "⠏ |         0 | task1                                                           \n",
            "\u001b[32m✔\u001b[39m |         3 | task1\n",
            "\n",
            "⠏ |         0 | task-2-subtask-1-tests                                          \n",
            "⠏ |         0 | task2a                                                          \n",
            "\u001b[31m✖\u001b[39m | \u001b[33m2\u001b[39m       0 | task2a\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\u001b[1m\u001b[33mFailure\u001b[39m (test-task-2-subtask-1-tests.R:6:3): Marking task2\u001b[22m\n",
            "exists(\"task2\") is not TRUE\n",
            "\n",
            "`actual`:   \u001b[32mFALSE\u001b[39m\n",
            "`expected`: \u001b[32mTRUE\u001b[39m \n",
            "task2() does not exist.\n",
            "\n",
            "\u001b[1m\u001b[33mError\u001b[39m (test-task-2-subtask-1-tests.R:7:3): Marking task2\u001b[22m\n",
            "Error in `eval_bare(expr, quo_get_env(quo))`: object 'task2' not found\n",
            "\u001b[1mBacktrace:\u001b[22m\n",
            "\u001b[90m 1. \u001b[39m\u001b[1mtestthat\u001b[22m::expect_is(task2, \"function\", info = \"task2 is not a function.\")\n",
            "\u001b[90m      at test-task-2-subtask-1-tests.R:7:2\u001b[39m\n",
            "\u001b[90m 2. \u001b[39mtestthat::quasi_label(enquo(object), label, arg = \"object\")\n",
            "\u001b[90m 3. \u001b[39m\u001b[1mrlang\u001b[22m::eval_bare(expr, quo_get_env(quo))\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\n",
            "⠏ |         0 | task-2-subtask-2-tests                                          \n",
            "⠏ |         0 | task2b                                                          \n",
            "\u001b[31m✖\u001b[39m | \u001b[33m1\u001b[39m       0 | task2b\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\u001b[1m\u001b[33mError\u001b[39m (test-task-2-subtask-2-tests.R:6:3): Mark even more on task2\u001b[22m\n",
            "Error in `task2(5:10)`: could not find function \"task2\"\n",
            "\u001b[1mBacktrace:\u001b[22m\n",
            "\u001b[90m 1. \u001b[39m\u001b[1mtestthat\u001b[22m::expect_is(task2(5:10), \"integer\", info = \"task2 don't return an integer.\")\n",
            "\u001b[90m      at test-task-2-subtask-2-tests.R:6:2\u001b[39m\n",
            "\u001b[90m 2. \u001b[39mtestthat::quasi_label(enquo(object), label, arg = \"object\")\n",
            "\u001b[90m 3. \u001b[39m\u001b[1mrlang\u001b[22m::eval_bare(expr, quo_get_env(quo))\n",
            "────────────────────────────────────────────────────────────────────────────────\n",
            "\n",
            "══ \u001b[1mResults\u001b[22m ═════════════════════════════════════════════════════════════════════\n",
            "\u001b[36mDuration: 0.8 s\u001b[39m\n",
            "\n",
            "[ \u001b[33mFAIL\u001b[39m 4 | WARN 0 | SKIP 0 | PASS 3 ]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "That did not work very well. This is how it looks when something goes wrong when the markmyassignment package is used. It looks daunting. The problem here is that `task2` is missing and also the mandatory variable `my_name` is missing.\n",
        "\n",
        "\n",
        "Below we fix the errors and show how it would look when everything works as it should. For a more indepth walk-trhough, see the full vignette [here](https://cran.r-project.org/web/packages/markmyassignment/vignettes/markmyassignment.html)."
      ],
      "metadata": {
        "id": "D4czGQuQQkw-"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Lets add task2\n",
        "task2 <- function(vector) vector[1] + vector[length(vector)]\n",
        "# And also add our name\n",
        "my_name <- \"Thomas Bayes\"\n",
        "# And lets check the results\n",
        "mark_my_assignment()"
      ],
      "metadata": {
        "id": "dTgcmjWMRXYR",
        "outputId": "b03af28d-e434-4377-87f7-db8d7a415ca1",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "\u001b[32m✔\u001b[39m | \u001b[33mF\u001b[39m \u001b[35mW\u001b[39m \u001b[34mS\u001b[39m \u001b[32m OK\u001b[39m | Context\n",
            "\n",
            "⠏ |         0 | mandatory-1                                                     \n",
            "⠏ |         0 | Mandatory tests                                                 \n",
            "\u001b[32m✔\u001b[39m |         1 | Mandatory tests\n",
            "\n",
            "⠏ |         0 | task-1-subtask-1-tests                                          \n",
            "⠏ |         0 | task1                                                           \n",
            "\u001b[32m✔\u001b[39m |         3 | task1\n",
            "\n",
            "⠏ |         0 | task-2-subtask-1-tests                                          \n",
            "⠏ |         0 | task2a                                                          \n",
            "\u001b[32m✔\u001b[39m |         3 | task2a\n",
            "\n",
            "⠏ |         0 | task-2-subtask-2-tests                                          \n",
            "⠏ |         0 | task2b                                                          \n",
            "\u001b[32m✔\u001b[39m |         4 | task2b\n",
            "\n",
            "══ \u001b[1mResults\u001b[22m ═════════════════════════════════════════════════════════════════════\n",
            "[ FAIL 0 | WARN 0 | SKIP 0 | \u001b[32mPASS\u001b[39m 11 ]\n",
            "Good work!"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Using RStan/Stan\n",
        "Below are two simple examples of running Stan for two simple models. First we just try to estimate the mean of a data where the standard deviation is known, then we take the classical `rats` example. "
      ],
      "metadata": {
        "id": "x9st_16UzZJw"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "### The unknown mean model"
      ],
      "metadata": {
        "id": "XkEY15LVORdu"
      }
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "lxBXXzeIx6aA"
      },
      "source": [
        "library(rstan)\n",
        "stan_model_code <- \"\n",
        "data {\n",
        "    int<lower=0> N; \n",
        "    real y[N];\n",
        "}\n",
        "\n",
        "parameters {\n",
        "    real mu; \n",
        "\n",
        "}\n",
        "\n",
        "model {\n",
        "    target += normal_lpdf(mu | 0, 10);\n",
        "    target += normal_lpdf(y | mu, 1);\n",
        "}\"\n",
        "\n",
        "\n",
        "# Create an example data set\n",
        "set.seed(4711)\n",
        "dat <- list(N = 20, y = rnorm(20));\n",
        "\n",
        "# Run rstan for the model\n",
        "fit <- stan(model_code = stan_model_code, model_name = \"example_model\",\n",
        "            data = dat, iter = 2000, chains = 4)\n"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "As the next step we check that we don't have any divergent transitions."
      ],
      "metadata": {
        "id": "kSqeqzk4L0bO"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# We check that no HMC draw did diverge\n",
        "check_divergences(fit)"
      ],
      "metadata": {
        "id": "t_352u9JLzw1",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "c75f6ec7-5fd7-43e9-be12-3114dfeb8e8f"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "0 of 4000 iterations ended with a divergence.\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "To look at the results we can simply use `print` in R."
      ],
      "metadata": {
        "id": "keQtRTzoGcse"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Print the results\n",
        "print(fit)"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "kZ67roHUGlC-",
        "outputId": "dd9c5107-67d5-4c0d-c541-94618c527a46"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Inference for Stan model: example_model.\n",
            "4 chains, each with iter=2000; warmup=1000; thin=1; \n",
            "post-warmup draws per chain=1000, total post-warmup draws=4000.\n",
            "\n",
            "       mean se_mean   sd   2.5%    25%    50%    75%  97.5% n_eff Rhat\n",
            "mu     0.12    0.01 0.22  -0.32  -0.03   0.11   0.27   0.55  1545    1\n",
            "lp__ -31.53    0.02 0.70 -33.49 -31.68 -31.26 -31.08 -31.04  1659    1\n",
            "\n",
            "Samples were drawn using NUTS(diag_e) at Sun Jun 26 09:29:55 2022.\n",
            "For each parameter, n_eff is a crude measure of effective sample size,\n",
            "and Rhat is the potential scale reduction factor on split chains (at \n",
            "convergence, Rhat=1).\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "If we don't have too many parameters, we can also visualize the results with `plot`."
      ],
      "metadata": {
        "id": "Cczs39_tGpN1"
      }
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Vmo1YXVcx91I",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 268
        },
        "outputId": "a90e1ae0-71ff-460b-b2c8-f4824993dc6a"
      },
      "source": [
        "plot(fit)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "ci_level: 0.8 (80% intervals)\n",
            "\n",
            "outer_level: 0.95 (95% intervals)\n",
            "\n"
          ]
        },
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "plot without title"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAFoCAIAAAD1h/aCAAAACXBIWXMAABJ0AAASdAHeZh94\nAAARB0lEQVR4nO3ca3BUZZ7H8acv6ZBOAqEJJKYDjCkkMorLRHGCDquj4jrgBSwpBYsAWhoW\nXV4oo4LDdRiYEC+wgIukoHDXWgpLhXJcvHBzlZWaWVQa14i4sSpGFEICSediSEh6X/QaEEhy\nfiHdpxu+n1fdycnpfx5Of3P6hiMUChkAUDjtHgBA/CEcAGSEA4CMcACQEQ4AMsIBQEY4AMgI\nBwBZzIXj5ZdfvvnmmwOBwOHDh+2eRRBH09bV1QWDQbunsCqOFtZcStO6e3COHlFTU1NRUXHy\n5MlTp07ZPYsgjqZta2uLo7cLx9HCmktp2pg74wAQ+wgHABnhACAjHABkhAOAjHAAkBEOADLC\nAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJAR\nDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CM\ncACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAjHABk\nhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAg\nIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcA\nGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgA\nyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIB\nQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEO\nADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxw\nAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSE\nA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAj\nHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ\n4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkLntHgAXlfqDZVW7\nPz75Q6VxOpMGXtZ/zOikQX67h0LPIxzoGdUf7C19etnxv+3/X9NcaU65jCPDuHMcnv5jRv+y\naE7vEVfZPSB6EuFAD/j6T6sO/GH526Z+iwnWmrb2rw8IuSe+/87vPvzr8DVLBj10v40TomcR\nDlyosufX7ftD0RJTVWpOnvWtSnNqjTn+adOPTz48252akjVxnC0Tosd1EY6NGze++eabxph5\n8+ZlZWWVlJSUlpYmJSXddNNNU6dOdblc77///pYtWyorKzMzMydPnnzjjTeGf7C4uPijjz4y\nxsydOzc/Pz/8xc8///zZZ581xuTl5S1cuDByvxWipu5/vjo4589FpvrcarTba378F3PC+/Dv\n+/39rxMz0qM5HiKki3B4PJ7whaNHj65Zs+b48ePGmKampq1bt7rd7uzs7NWrV4c3qKioWL58\n+fPPPz9kyJCIThx3av47cOAf59o9xWmtra3GGJfL1SN7a/zm2x0ttQHT1PlmO03DrXXVzr+7\nPSn7Mmn/zc3N5T8dhLGv72NTzPTBdk8RDV2Eo/3weuONN/x+f0FBwa5duw4cOGCMee+99xIS\nEsaNG5eTk7Np06aqqqpQKLRt27ZZs2ZFfOq4cqquofaTz+2eIoK2mXqLm119tKr5aJW6/x/1\nkezS+0St3SNEidXnOFpaWubPn5+YmHj99ddPmTKltbW1rq4uPz+/sLDQGOP1eouKiowx33zz\nTTeGaG5uPnbsWPhyQ0NDN/bQiZqamvCJUkQdPnw4/Jf8XA0//BDpW7dRkwl9bZqtbHmgq7OS\ni0B1dXX37gI9wuv1ZmZmRue2rIZj5MiRiYmJxpiUlJTs7Ozy8nJjzOjRo8PfHTZsWPhC9+6i\nhw4dmjZtWvtVj8fz/fffOxyOju6KkjVr1rzwwgsXvp9uu8b0WmL62zhARJ0wrSFrWwZNW4sJ\nJRhHZAey1dKlS3ctfdauW7/llltKSkqsb19dXd3JXcztdg8aNKjD71q8jfT0089ppaSkhC/0\n79//rK90767u8/nuvffe8OX9+/cHAoGsrKwBAwYMHtwDDxdvvfXW+npL59IXor6+vn0RzpL+\n/XHz9l8jPYBdelkOQYJxuC/qahhjRo0aNWR4jl23Pnz48Jwc4dZdLle372JWw+F2n97S4XCc\n+8WOhEKn/yC1tLScd5usrKy5c///6cOioqJPPvnE4lRWjB07duzYsT24w/MqLy/v6N+gatfH\ne99+INID2CXNuPoY55nv3ejIIJNwkWfDmPvuu2/EE4V2TxENkfqsSvuzqmf+tT9y5EiEbg52\ncRhzg/Fa2fI3JinSwyBqIvUGMJ/PF74QCATGjBljjGlra9u+fXuEbi6W9Rnxy/zt/273FKc1\nNDSEQqGOHlipjv5lx8R/LtltGppMZ8919DOucSb1yj/9Pu36X2n7P3o0IyPjwmaMnto+lhp6\nEYhUONrfzfHhhx/6fL4hQ4bs2rUrGAxG6OZiWYIvrf9tv7F7itM8tbWhUCgtLa1H9tbvpvwT\nH+97Yl9zkanq6PmtRON4xqT/4p5/uGLuP6n7bywv798TT3VFR2N5ud0jREmkHqqMGjWq/ZWh\nrVu3Pvfcc1988cWjjz4a/kpbW9cPiREXnAnukVtKbrl86ELTP92c501l2SahyAy4Pi/vV/+6\nIvrjIUIidcbhcrmWLFmyfv36QCAQCoWGDh1aUFCQlZUV/m5T08X/kv6lo1f2ZaP/9pfUaU8M\n+4+d/2ka95kfj5pWlzGXGfevTdINJukXUycOX7PElXypnMZfChxnvuoRC4qKitatW7dp06aM\njIweeTk2Ojp5VSXW1PboQ5UzVe3++Lt/e7Nq13+d/P6ocTqTBvn7jxk96OEH+uRd3e19xtHC\nmktpWj4dix6T/tsb0n97g91TIBr4rwMByAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADI\nCAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFA\nRjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4A\nMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAA\nkBEOADLCAUBGOADICAcAGeEAICMcAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQD\ngIxwAJARDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUBGOADICAcAGeEAICMc\nAGSEA4CMcACQEQ4AMsIBQEY4AMgIBwAZ4QAgIxwAZIQDgIxwAJARDgAywgFARjgAyAgHABnh\nACAjHABkbrsHOL+33norJSWlb9++dg9i1YkTJ+Jl2qamplAolJSUZPcglsTRwpqLa9revXtP\nmjSpo+/GaDg2bdpk9wjAJW3gwIGdhCPmHqrcc8891113ndvtTk5OdrtjtGtnSUxMTE5Odjgc\ndg9iidfr9Xq9dk9hidPpTE5OTkxMtHsQSxISEuLooO3Vq9eFHLQxF44rr7wyKyvL6XR6PB6X\ny2X3OJa43W6Px+N0xtxinldCQoLH47F7CkscDofH44mXu2I8HrTdDocjFAr17EAXrrGx8Z13\n3lm2bNnjjz8+YcIEu8fpWnFx8bvvvrthw4bBgwfbPUvXHnjggVAotHnzZrsH6VpFRcW0adNu\nv/32p59+2u5ZurZ169ZVq1Y988wzY8aMsXuWri1YsGDPnj2vvfZav379zruB0+lMTU3t6Mdj\nseXt59JJSUl9+vSxe5yuhf+Ap6amxsW0TqczFArFxag1NTXGGI/HExfThp9v9nq9cTFtQkKC\nMaZ3797dmzYWw2GM8Xq9fr8/JSXF7kEsSUtL8/v94X+J2JeZmRmDp5nn5Xa7/X5/vLxOkZKS\n4vf74+X5I5/P5/f7u/34OhYfqgCIcfHxfB6AmEI4AMgIBwBZrDw5Wl9fv27dugMHDrS0tOTm\n5s6YMWPAgAEdbbxz586VK1fOnTs3Pz8/mkO2szLt8ePHN2zYEAgEmpubc3Jypk+fPnTo0Fgb\nUlp226e1fUmtj9ouLo5VY8y2bdu2bNlSXV3t9/sLCgpGjhzZ+W5j5YxjxYoVlZWVCxYsKC4u\n9nq9ixcvbmtrO++WNTU1r7zyir1vYbIy7ZIlS6qqqhYtWrRixYr09PTFixc3NTXF2pDWlz0W\nprV9Sa2PGhYvx+rOnTs3b95cWFi4du3a2267raSkpLGxsYv9hmLAsWPH7r777rKysvDVurq6\n8ePH79+//7wbL1u2bP369VOmTNm7d28UZzzNyrTBYHDp0qXffvtt+GplZeVdd9116NChmBpS\nWvaIiosltT5qu7g4VkOh0COPPLJz505pzzFxxvH1118nJCRcfvnl4aspKSnZ2dlfffXVuVvu\n3bu3rKxs8uTJ0R3wZ6xMm5qaOmfOnIEDB4avVldXO53O9PT0mBrS+rJHWlwsqfVRw+LlWK2u\nrj5y5IgxZtasWRMnTpw9e/bBgwe73HNMhCMYDKampp75tvk+ffrU1taetVl9ff3atWsfe+yx\nXr16RXfAn7E4bbu6urpVq1aNHz8+mm9ksjKk+otETlwsadjFd6xWV1cbY3bs2PHUU09t2LAh\nNzd30aJFXR4G9oRjz54943/y5ZdfGmOsfNhm/fr1eXl5I0aMiPyAP9O9acO+++672bNnX331\n1VOnTo3kjOdhZcjY+URvXCxpWCwfq+eyuLD3339/dnZ2amrqQw895HA49u3b1/n29ryqkpeX\nt3LlyvDlzMzMYDAYDAZDoVD7L1lbW3vWH5P9+/d/+umnq1evjvas3Zo2LBAILF++fNKkSXfe\neWf0xjXGGJOWltblkFa2iQ7rk9i4pGFWRrXxWD2LlWl9Pp8xJjk5OXzV5XL5fL4TJ050vmd7\nwuH1es/8IOkVV1zR0tJSVlY2ZMgQY0wwGKyoqBg2bNiZP7J9+/aGhoYZM2aEr9bX17/44osj\nRoyYM2dODE5rjCktLS0qKnryySevvfbaSE94LitDWvxFYmRaY/eShsX4sdqNaX0+X9++fQ8e\nPBjeprm5+dixYxkZGZ3v2bVw4cKIjW1VUlJSeXn57t27c3NzGxsbX3rppeTk5AcffNDhcGzf\nvr20tDQ3N/eaa6753Rk++OCD6dOnT5gwIfr/y4uVaZubm+fPn3/HHXfk5eU1/sTpdEbtv5aw\nMmQn20RnSGla25fU+qjxdaw6HI7W1tbXX389JyfH7XZv3LixsrKysLCw84WNlQ+5NTY2rlu3\n7rPPPmttbb3qqqtmzJgRPqEqLi4OBoN//OMfz9q+oKBg5syZdr2ppstpA4HAvHnzzvqpwsLC\ncePGxc6QnWwTfXGxpBZHPWv7GD9WjTFtbW2vvvrqjh076uvrc3NzZ86c2f7qVUdiJRwA4khM\nvBwLIL4QDgAywgFARjgAyAgHABnhACAjHABkhAOAjHAAkBEOADLCAUD2f7Y8WvxsmUhMAAAA\nAElFTkSuQmCC"
          },
          "metadata": {
            "image/png": {
              "width": 180,
              "height": 180
            }
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "### The rats model\n",
        "Below we run a slightly more complicated model. We start out by setting up the data."
      ],
      "metadata": {
        "id": "0wW8WTNvJuhz"
      }
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "0ocs9zTCx-Sw",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "5f753993-cbd3-47ff-83b6-38cddfe9dccb"
      },
      "source": [
        "# Create the rats dataset\n",
        "y <- as.matrix(read.table('https://raw.github.com/wiki/stan-dev/rstan/rats.txt', header = TRUE))\n",
        "x <- as.array(c(8, 15, 22, 29, 36))\n",
        "xbar <- mean(x)\n",
        "N <- nrow(y)\n",
        "T <- ncol(y)\n",
        "dat <- list(N=N, T=T, y=y, x=x, xbar=xbar)\n",
        "print(dat)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "$N\n",
            "[1] 30\n",
            "\n",
            "$T\n",
            "[1] 5\n",
            "\n",
            "$y\n",
            "      day8 day15 day22 day29 day36\n",
            " [1,]  151   199   246   283   320\n",
            " [2,]  145   199   249   293   354\n",
            " [3,]  147   214   263   312   328\n",
            " [4,]  155   200   237   272   297\n",
            " [5,]  135   188   230   280   323\n",
            " [6,]  159   210   252   298   331\n",
            " [7,]  141   189   231   275   305\n",
            " [8,]  159   201   248   297   338\n",
            " [9,]  177   236   285   350   376\n",
            "[10,]  134   182   220   260   296\n",
            "[11,]  160   208   261   313   352\n",
            "[12,]  143   188   220   273   314\n",
            "[13,]  154   200   244   289   325\n",
            "[14,]  171   221   270   326   358\n",
            "[15,]  163   216   242   281   312\n",
            "[16,]  160   207   248   288   324\n",
            "[17,]  142   187   234   280   316\n",
            "[18,]  156   203   243   283   317\n",
            "[19,]  157   212   259   307   336\n",
            "[20,]  152   203   246   286   321\n",
            "[21,]  154   205   253   298   334\n",
            "[22,]  139   190   225   267   302\n",
            "[23,]  146   191   229   272   302\n",
            "[24,]  157   211   250   285   323\n",
            "[25,]  132   185   237   286   331\n",
            "[26,]  160   207   257   303   345\n",
            "[27,]  169   216   261   295   333\n",
            "[28,]  157   205   248   289   316\n",
            "[29,]  137   180   219   258   291\n",
            "[30,]  153   200   244   286   324\n",
            "\n",
            "$x\n",
            "[1]  8 15 22 29 36\n",
            "\n",
            "$xbar\n",
            "[1] 22\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "We now load the rats Stan model code. Note that `//` indicate comments"
      ],
      "metadata": {
        "id": "oyo_GFKWKMRj"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# rats_stan_code <- readLines('https://raw.githubusercontent.com/stan-dev/example-models/master/bugs_examples/vol1/rats/rats.stan')\n",
        "rats_stan_code <- \"\n",
        "// http://www.mrc-bsu.cam.ac.uk/bugs/winbugs/Vol1.pdf\n",
        "// Page 3: Rats\n",
        "data {\n",
        "  int<lower=0> N;\n",
        "  int<lower=0> T;\n",
        "  real x[T];\n",
        "  real y[N, T];\n",
        "  real xbar;\n",
        "}\n",
        "parameters {\n",
        "  real alpha[N];\n",
        "  real beta[N];\n",
        "  \n",
        "  real mu_alpha;\n",
        "  real mu_beta; // beta.c in original bugs model\n",
        "  \n",
        "  real<lower=0> sigmasq_y;\n",
        "  real<lower=0> sigmasq_alpha;\n",
        "  real<lower=0> sigmasq_beta;\n",
        "}\n",
        "transformed parameters {\n",
        "  real<lower=0> sigma_y; // sigma in original bugs model\n",
        "  real<lower=0> sigma_alpha;\n",
        "  real<lower=0> sigma_beta;\n",
        "  \n",
        "  sigma_y = sqrt(sigmasq_y);\n",
        "  sigma_alpha = sqrt(sigmasq_alpha);\n",
        "  sigma_beta = sqrt(sigmasq_beta);\n",
        "}\n",
        "model {\n",
        "  mu_alpha ~ normal(0, 100);\n",
        "  mu_beta ~ normal(0, 100);\n",
        "  sigmasq_y ~ inv_gamma(0.001, 0.001);\n",
        "  sigmasq_alpha ~ inv_gamma(0.001, 0.001);\n",
        "  sigmasq_beta ~ inv_gamma(0.001, 0.001);\n",
        "  alpha ~ normal(mu_alpha, sigma_alpha); // vectorized\n",
        "  beta ~ normal(mu_beta, sigma_beta); // vectorized\n",
        "  for (n in 1 : N) {\n",
        "    for (t in 1 : T) {\n",
        "      y[n, t] ~ normal(alpha[n] + beta[n] * (x[t] - xbar), sigma_y);\n",
        "    }\n",
        "  }\n",
        "}\n",
        "generated quantities {\n",
        "  real alpha0;\n",
        "  alpha0 = mu_alpha - xbar * mu_beta;\n",
        "}\"\n",
        "cat(rats_stan_code, sep = \"\\n\")\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "sOUi6esxH8C6",
        "outputId": "f2a098be-cdb1-4078-c416-74b11a7b411a"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "\n",
            "// http://www.mrc-bsu.cam.ac.uk/bugs/winbugs/Vol1.pdf\n",
            "// Page 3: Rats\n",
            "data {\n",
            "  int<lower=0> N;\n",
            "  int<lower=0> T;\n",
            "  real x[T];\n",
            "  real y[N, T];\n",
            "  real xbar;\n",
            "}\n",
            "parameters {\n",
            "  real alpha[N];\n",
            "  real beta[N];\n",
            "  \n",
            "  real mu_alpha;\n",
            "  real mu_beta; // beta.c in original bugs model\n",
            "  \n",
            "  real<lower=0> sigmasq_y;\n",
            "  real<lower=0> sigmasq_alpha;\n",
            "  real<lower=0> sigmasq_beta;\n",
            "}\n",
            "transformed parameters {\n",
            "  real<lower=0> sigma_y; // sigma in original bugs model\n",
            "  real<lower=0> sigma_alpha;\n",
            "  real<lower=0> sigma_beta;\n",
            "  \n",
            "  sigma_y = sqrt(sigmasq_y);\n",
            "  sigma_alpha = sqrt(sigmasq_alpha);\n",
            "  sigma_beta = sqrt(sigmasq_beta);\n",
            "}\n",
            "model {\n",
            "  mu_alpha ~ normal(0, 100);\n",
            "  mu_beta ~ normal(0, 100);\n",
            "  sigmasq_y ~ inv_gamma(0.001, 0.001);\n",
            "  sigmasq_alpha ~ inv_gamma(0.001, 0.001);\n",
            "  sigmasq_beta ~ inv_gamma(0.001, 0.001);\n",
            "  alpha ~ normal(mu_alpha, sigma_alpha); // vectorized\n",
            "  beta ~ normal(mu_beta, sigma_beta); // vectorized\n",
            "  for (n in 1 : N) {\n",
            "    for (t in 1 : T) {\n",
            "      y[n, t] ~ normal(alpha[n] + beta[n] * (x[t] - xbar), sigma_y);\n",
            "    }\n",
            "  }\n",
            "}\n",
            "generated quantities {\n",
            "  real alpha0;\n",
            "  alpha0 = mu_alpha - xbar * mu_beta;\n",
            "}\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        ""
      ],
      "metadata": {
        "id": "njlBOHrgIgK_"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Run rstan for the rats model\n",
        "rats_fit <- stan(model_code = rats_stan_code, model_name = \"rats_model\",\n",
        "                 data = dat, iter = 2000, chains = 4)"
      ],
      "metadata": {
        "id": "KN-xMUzxIg3M"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# We check that no HMC draw did diverge\n",
        "check_divergences(rats_fit)"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "mNPF44h0LlHI",
        "outputId": "20ef16f0-3615-4edc-da52-d0094d0f5b4d"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "0 of 4000 iterations ended with a divergence.\n",
            "\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "# Print the results\n",
        "# Here we just print a few parameters\n",
        "print(rats_fit, pars = c(\"beta[1]\", \"mu_beta\"))\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "kBdzx96tKToy",
        "outputId": "e4c1a0ca-4926-492c-dde1-e197371ec6d0"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Inference for Stan model: rats_model.\n",
            "4 chains, each with iter=2000; warmup=1000; thin=1; \n",
            "post-warmup draws per chain=1000, total post-warmup draws=4000.\n",
            "\n",
            "        mean se_mean   sd 2.5%  25%  50%  75% 97.5% n_eff Rhat\n",
            "beta[1] 6.06       0 0.25 5.56 5.90 6.06 6.23  6.54  5758    1\n",
            "mu_beta 6.19       0 0.11 5.98 6.11 6.19 6.26  6.40  4418    1\n",
            "\n",
            "Samples were drawn using NUTS(diag_e) at Sun Jun 26 09:30:51 2022.\n",
            "For each parameter, n_eff is a crude measure of effective sample size,\n",
            "and Rhat is the potential scale reduction factor on split chains (at \n",
            "convergence, Rhat=1).\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "# Print the results\n",
        "plot(rats_fit, pars = c(\"beta[1]\", \"mu_beta\"))"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "_oREzqELKmPb",
        "outputId": "fe8990b8-9bb7-42ff-ea65-f63bb1459c41"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "ci_level: 0.8 (80% intervals)\n",
            "\n",
            "outer_level: 0.95 (95% intervals)\n",
            "\n"
          ]
        },
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "plot without title"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAFoCAIAAAD1h/aCAAAACXBIWXMAABJ0AAASdAHeZh94\nAAAd9UlEQVR4nO3ce1zUdb748fdcgGGGu4gmt5VQrG1zj6mbdlmtR6dO5u6mbbi2GrWa1Vab\nqZiobRt2OqbJrpm27pbmj8JsTda00khq+VmetEUCL2sqISA3FQGBgWGY88fXZgkF5jMCw+Tr\n+U8zw2eG9xjzYr6XQedwOAQAVOg9PQAA70M4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxw\nuOnPf/7zuHHj8vLyPD3Id9jt9vLyck9P4abS0lJPj+Cmy3BywuGms2fPFhcXNzU1eXqQ9lpa\nWjw9gpuYvPe5PTnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgj\nHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCM\ncABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4Ay\nwgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDK\nCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAo\nIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCg\njHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOA\nMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4A\nyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQJnR0wMArrLXN1jL\nKg3+Jt+IcL0PP7qe1CP/+m+//fabb74pIvfdd19iYmJPfIsuzZ0798iRI86rJpNp06ZNXd5r\n8eLFeXl5bW/ZunVr9w8HFfaGxm9Wbyh9K7Mm94B2i8HfFH7rDT94dHrEf4337GyXLS/IdnZ2\ndlpa2k033TRv3rxuecDm5uaMjIwtW7a0traKyN/+9jdfX99ueWR0u+rPv9z3y0eOlhZvl3P7\nxVoldj/RRTb6jN1Wdce2rNgJt/5H+kqfkCBPj3nZ8YJw7Nmz51Luvm7dOh8fH51Op109cODA\nqlWrSktLL7p4wYIFLS0tIrJ48eLCwsJL+b64dKc+3v2/dyW9aa3aJDX2b29sFsdhaTosTVuk\nbv72DxpunHTD7i0+wYGeHPTy09d3jjY3N+fm5l7KIwQGBgYFBQUGBopIVlZWSkpKaWlpZGSk\nj4/PhYvNZnNQUFBQUJDBYLiUb4pLZy0t33fvI3+1VmS0qUZb1WJfLFV7D+TvT3qqt4e77PXs\nOw69Xl9aWrphw4aCggKbzRYfHz9t2rSrrrrKucDhcOzatSsrK6uwsLClpSUiImLcuHGTJk0y\nGo3y3T0OOTk5OTk5MTExq1atEpGGhobNmzfv2bOnoqLC4XAMHDjwxhtvnDx5cucbHeXl5QaD\nYfLkyffee++0adNsNltPPvs+p2rHpzX7D3p6ig5VV1fbQkOdV8s2f7D3THmm1HVyF5s4XpTT\ncZkf+j70tOXK2J6f8eLaTe6KyF/9zD8msofm6QU9G46zZ88mJyfX1Z3/f19QULBo0aKlS5fG\nx8eLiMPhWL58eU5OzvlRjMaSkpL09PS8vLznnnuuk9/5TU1NKSkpx48fd95SXFyckZFRUFCw\nZMkS51bJhYYMGbJy5cqoqKjueXrepixzR9Gr6Z6eojPl3726UWq7vMspsWdJve9f3uqhkVxU\n3vWS7wj9yY+9Ohw9u6myY8eO4cOHr1ixIjU1NTY2VkRsNtv69eu1r2ZlZWnViI+PX79+/Tvv\nvDN16lQRyc/Pf//990Vk4cKFTz/9tLZ4zJgxGRkZy5YtE5HMzEytGqNHj964cePrr78eGRmp\n3XHv3r2dzPOTn/zksq2G16kR+yFpcmXlHmns6WHQTs++4wgKCpo9e7a2N2Hu3LmPP/64iOTn\n59fX11sslu3bt2vLkpKSwsLCRCQxMXHXrl3l5eUffvjhxIkTTSaTyWQ6P6jRaLFYnJdvueUW\nEbnrrrvMZrPZbB4/fnx6erqIHDx4cPTo0T30dI4dO/bss89ql8vLy41Go81ma25udu/RsrKy\nFi9e3G3DiYiIw+FoaWm56O4bEZlwom5U936/nlQmdodrK0vF+zY5H3rooW8Ce/ZYXmJi4pNP\nPtn5mpaWlo5+gHU6XUc/SNLT4Rg1apTze8fGxoaHh586dcrhcFRUVAwaNMh52ELbctFmHTZs\nWHl5eXFxcU1NTXBw8EUfdvLkye1u0bojIo2NPfjLp7Gx8dChQ86rer3+3Llz1dXV7j1acXHx\nP//5z24azSU/kVCRgN78jpeiRVzshrSc/69DpMOt1L7m66+/znft/ZTbRo4c2eUPZ319fUdr\nDAZDeHh4R3fs2XAMGDCg7dXQ0NBTp06JyJkzZwICAhyO8z8ZU6ZMufC+J0+e7CgcNptt586d\n2dnZZWVlzh0oGudj9oRrrrlm37592uWlS5euXbs2NDS03XN03axZs2bNmtV904mI2O32srKy\njjbHvnokpY/v42irv7h6YCvi/EqvqYaIZGdn9xs3xtNTiNVqde8HuFfP43Du79TOlXCKiIi4\ncHG7NW298MIL2gvYbDaPGDHCZDJVVlYePXq0W4eFZzkGiDFSjKXS4Y+B0wjx74WB0FbPhqOq\nqqrt1TNnzmgXQkJCQkJC9Hq9du7m6tWrXT9389ixY1o1LBbLmjVrQkJCRGTbtm2Eo0vBI64Z\n9MsJnp6iQw0NDWaz2Xn1TM4Xd5fXr5Iznd/LJLo7xBI6+sf+sR47SNFuclf49u/XQ8P0jp4N\nx759+2bMmKG90SgpKamsrBQRg8EQGxvr6+sbFxenvdoPHz587bXXanfRDpcEBQX169ev7YFV\n5y6c8vLzR76io6O1aoiI8ySxHt1U8XaxM6fGzpzq6Sk6VFRUpB1609T8s6Duuv/6TBr+KdZO\n7jVTQmOvveaGz7boPHfOXrvJLwc9G47KysoVK1YkJibW1dWtXr1au3HkyJH+/v4icuedd65c\nuVJEVq9e/dhjj8XExOTn56elpTU3N0dHR2snejmPpBQUFBw4cMBisfTv31+7paioqLS0NCQk\nZNu2bQUFBTqdzuFwFBYW2u12vf7ih5m/+eYb57se7c2OiOzfv1873ywqKuqiG03wiOAR1/z4\nj39IfvKZZXL6y4u1wyDyoITcGXLFyHfWeLAal6ceCYdz98Ttt9+enZ3tPMVLRMxmc1JSknb5\n1ltv/fLLL3fv3n3y5MmUlBTnmoCAgOTkZO3tRkxMjMViqa+vr6+vX7BggclkysjIiIuLO378\neGNj4yOPPKItmzNnzqZNm06cOHHkyJFJkyalpaVddLB33333k08+aXfjkiVLtAsPPPDA3Xff\n3R3/AOgecb97UGfQ/2FO6sfNZ7fKuWNy/l2nr+iuE9MUCf7hDwaPyvyrZWicZ+e8DPVIOJyb\nFQkJCXfccccbb7xx6NAhvV5/9dVXJyUlaSdriYhOp0tOTt65c2dWVtaJEyeam5v79es3YsSI\ne+65x/mb39/ff968ea+99lpZWZm/v39CQoLBYFi0aNG6devy8/OtVmtsbGxiYuLIkSNDQ0Nf\nfvnlysrKyMjIgACvOeiIzg1+LCl8/NiIxctv3f7x2eam02I3igwQY0D/8B88/Osr580yBvL/\n2gN039edAs6/x+Hep+afeuopbf9LR3+PQzscm5GR0XPnm7mh88OxfVznewpsNXVn/zfXWlap\n9/O1xMUEX/ejvrN54r37ONye3As+Vg+IiE9wYP//vNnTU+C873846urqtL/HoX2yvnMNDQ3a\nDhq7/aKf5AYgcjmE44EHHhCX/3TgCy+80O5PBwK4UF//Qz4A+qDv7TuO5cuXu3Gv1NTUbp8E\n+P7hHQcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBG\nOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ\n4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBl\nhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACU\nEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQ\nRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFA\nGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKDN6eoCuvf3222++\n+aaI3HfffYmJiZ4eB16vtam5auc/qj7KaSwqcbS2mgYN6Dd+7IAJtxgDAzw9mtfwgnB0r+zs\n7LS0tJtuumnevHmengUeULpx66H5L5SfKP5SGsulpVUkXAzXrf1/g8L6DV38u7jfPSg6nadn\n9AKXXTj27Nnj6RHgMQfnLjnw0p83yNkPpN4uDuftOqm+8czZGbOfOZPzxYiMVXpfHw8O6RUu\nr30czc3Nubm5np4CnnF06Zrcl16dLxXb5FzbaoiIQyRHGuZIxb5338t/JMVTE3qRLt5xrF+/\n/t133xWRxYsXDxo06C9/+cvBgwf9/f1/+tOf3n///QaDYefOnVu2bKmsrBw4cODUqVNvuOEG\n7Y7Lli3LyckRkZSUlOuvv167MT8/f+HChSIyYsSIZ599VnVWvV5fWlq6YcOGgoICm80WHx8/\nbdq0q666yrnA4XDs2rUrKyursLCwpaUlIiJi3LhxkyZNMhqN2lPIy8vTVubk5OTk5MTExKxa\ntUpEGhoaNm/evGfPnoqKCofDMXDgwBtvvHHy5Mm+vr6qQ6Jvqv+68F/PLH9RTheKraM1p8S+\nRKr++PrGQYkT+//nzb05ntfpIhzOV05FRcUrr7xy5swZEbFarZmZmUajMSoqSnvhiUhxcfGL\nL7740ksvxcfH99CsZ8+eTU5Orqur064WFBQsWrRo6dKl2nd0OBzLly/XaiUiRqOxpKQkPT09\nLy/vueeeMxgMHT1sU1NTSkrK8ePHnbcUFxdnZGQUFBQsWbJExxZvn7HdFN/a1Oz23T+Txq/E\n2vmaUml5T+r8bv+16oN/9e0Fg7/pzoYj6tN5mS42VZyvt82bN0dGRj755JPXXnutdsuOHTs2\nbNgwYcKExx9/PDw8XEQcDsf777/fc7Pu2LFj+PDhK1asSE1NjY2NFRGbzbZ+/Xrtq1lZWVo1\n4uPj169f/84770ydOlVE8vPztakWLlz49NNPa4vHjBmTkZGxbNkyEcnMzNSqMXr06I0bN77+\n+uuRkZHaHffu3dtzTwe9LEvqu3HZZc7VnaM2m+2ZZ57x8/MbPXr0tGnT7HZ7XV3d9ddfP2vW\nLBExm81Lly4Vkba/t7tdUFDQ7NmzfXx8RGTu3LmPP/64iOTn59fX11sslu3bt2vLkpKSwsLC\nRCQxMXHXrl3l5eUffvjhxIkTTSaTyWTS1hiNRovF4rx8yy23iMhdd91lNpvNZvP48ePT09NF\n5ODBg6NHj3YOUFtb+8UXX2iXS0pK9PrLaw9Rtzty5Ihz41FEqqqq+vfv38l6U2vrpXy7w9Lk\nyrKT0lIjrcHu7v6z2+3vvPOO6+tvvvnmAQMGuPe9PMjVcIwaNcrPz09EAgICoqKiioqKROSm\nm27Svurc0aBty/SQUaNGadUQkdjY2PDw8FOnTjkcjoqKikGDBhUWFmpfcm4r6XS6YcOGlZeX\nFxcX19TUBAcHX/RhJ0+e3O4WrTsi0tjY2Pb2EydOON+ziIjRaKysrCwuLr7kZ9ZtWltbq6ur\nHQ5H10v7gPT09NTUVNfXb5YoH3Fzy7FVpF5c7U6t2N0OR3Nz87333uv6+vT09Jtv9tj+lNOn\nT3f0+89oNF5xxRUd3dHVcGgbI5qAgPPnyTh/PzhvsdvtLj6gG9qFOTQ09NSpUyJy5syZgIAA\n56tlypQpF9735MmTHYXDZrPt3LkzOzu7rKzMuQNF0+4VOGDAAO1tjoh8+umnOTk5ERER0dHR\n7j6h7me32w0GQ1RUlKcHccnEiROd7/tEpLq6OjQ0tJP1xkWvSoubP2B6EX/RN7jWjoBLONro\n4+PzP6n/4/r6sWPHevBHqLW11b3v7mo4tAMTGuf+wrY3dqTta89m63CHthuc+19aWlra3h4R\nEXHh4nZr2nrhhRf27dsnImazecSIESaTqbKy8ujRoxeu7N+///33369dLi8v//TTT90eHiIy\natSoUaNGOa8WFRVpu646sv33f2l1NxwiMkR887raOSoiEWIMlQ53pXfJaDTOnz/f7bt7i546\nAcz5qj537pzzxvLy8kt5zKqqqrZXnZtFISEhISEher2+tbVVRFavXu36YdRjx45p1bBYLGvW\nrAkJCRGRbdu2XTQc8Gq3iNmVcIwXcy8M4+16KhzO3QR5eXm33XabiLS2tn700UeX8pj79u2b\nMWOGlqSSkpLKykoRMRgMsbGxvr6+cXFx2qv98OHDzkM/2s7aoKCgfv36tT2w2tx8/qies2XR\n0dFaNUTEeZKYt+wsuEyMP5Qt7v4fOXf4qH3C/dvk3NfS2QHdcDHcLYHXrvlvpfM4SktLtSNx\nInKZnLHeU+Fw7qH8xz/+ERYWFh8fv2vXrtra2kt5zMrKyhUrViQmJtbV1a1evVq7ceTIkf7+\n/iJy5513rly5UkRWr1792GOPxcTE5Ofnp6WlNTc3R0dHa+ebOLeoCwoKDhw4YLFYnLtpioqK\nSktLQ0JCtm3bVlBQoNPpHA5HYWGh3W7X6/WczdEXmAe7vy/AHBeTkPzoghdXLZDKCrn4dmug\n6FMk/MpJE2IfVjuPw9fgMMfGuD2bN+qpcIwZM2bgwIHa7/PMzEwR8fPzmzt37vPPPy8irSrH\n1Zy7J26//fbs7GznKV4iYjabk5KStMu33nrrl19+uXv37pMnT6ak/Pus4YCAgOTkZO2VHxMT\nY7FY6uvr6+vrFyxYYDKZMjIy4uLijh8/3tjY+Mgjj2jL5syZs2nTphMnThw5cmTSpElpaWlx\ncXGX9u8Bzxv238nWkxUvpf/tz1L9/6Wh3VuX4WJ6VEKvGXfTj99Y4Zn5vEoP7uNYsmTJa6+9\nlpeX53A4hg4dOn369EGDBmlftVq73tR0cm5WJCQk3HHHHW+88cahQ4f0ev3VV1+dlJTkfIuo\n0+mSk5N37tyZlZV14sSJ5ubmfv36jRgx4p577nHuLvX39583b95rr71WVlbm7++fkJBgMBgW\nLVq0bt26/Px8q9UaGxubmJg4cuTI0NDQl19+ubKyMjIy0nnMCF5NZzD8x4a0kNHD+z2b9s2Z\nU19I40lpsYtjgBivE1OCf1Dc7BlDn31K73PZffLTDTo2492zdOnStWvXZmRktD1DzOPsdntZ\nWZm3HI5tp8ujKt3Fdra2bPP7VTv/0Vhc5rDZTFFXhI8fe8U9d5oGuXkiVq9N3u3cnpy44rLj\nExIU85spMb+5yPk+cJHHwlFYWPjiiy92uWz27NlDhw7thXkAuM5j4Rg8ePCaNWs89d0BXAo+\npgVAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDK\nCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAo\nIxwAlBEOAMoIBwBlhAOAMsIBQBnhAKCMcABQRjgAKCMcAJQRDgDKCAcAZYQDgDLCAUAZ4QCg\njHAAUEY4ACgjHACUEQ4AyggHAGWEA4AywgFAGeEAoIxwAFBGOAAoIxwAlBEOAMoIBwBlhAOA\nMsIBQBnhAKCMcABQZvT0AN5t69atubm5np7i31pbW+vq6oKDgz09iDuqq6tDQ0M9PYU7vpeT\nBwUF/epXv+rojoTjkmRkZHh6BKBHREdHdxIONlXc9POf/3zkyJGenqI9f39/i8Xi6Snc4efn\nZ7FYDAaDpwdR5uPjY7FYjEbv+x1sMBgsFouvr68b9yUcbho2bNigQYM8PUV7RqPR19dXp9N5\nehBlBoPBvZ9gj9Pr9b6+vt6YPJ1O5/bkOofD0e0DXSYaGhpsNpunp/iOp556Ki8v74MPPvC6\nF+HKlSv//ve/v/rqq0OGDPH0LGoyMzNffvnlp59++rbbbvP0LGr2798/Z86cqVOn/uY3v7nw\nq3q9PjAwsKP7et/7q77DbDZ7eoT2tDfMwcHBXhcObeCAgACv27Pr7+8vImaz2esm17Zq/fz8\n3JicTZXvlfDw8MjISE9P4Y6QkJDIyEiv652IBAQEREZG9sHfIl3y8/OLjIzs5G1FJ9hUAaCM\ndxwAlBEOAMoIBwBlHFXxSk888cQ333zjvGoymTZt2tR2QX5+/sKFC9vda9asWRMmTOiF8Tr3\n/vvvb9my5fTp05GRkdOnTx81alS7BefOnVu7du1XX31ls9kSEhIefvjhiIgIj4zaTpeTFxcX\nv/HGG4cOHXI4HIMHD542bdqwYcM8MuqFuhzexTXnOeCFHnjggffee6/qW6dPn263oLm5uaqN\nAwcO/PKXvzxx4oRHpm0rKytr+vTpe/furaioyMzMnDlzZn19fbs1qamp8+fPP378eGlp6bJl\ny37729/a7XaPTNtWl5PbbLb7779/xYoVJSUlJ0+eTEtLmzJlSkNDg6cGbsuVf3ZX1jgRDq90\nzz337N271/X1ixYteuutt3puHtfNnDnz448/7mRBVVXVz372s2PHjmlX6+rqfvGLX+zfv79X\nputMl5OfPXv23XffdZaipKRk4sSJzifiWV0O7+IaJzZVvI/NZmtqavr888/T09Pr6uri4+On\nT5/eyekbOTk5ZWVlv//973tzyIs6ffp0eXm5iDzxxBNlZWWxsbEzZsxo92b+66+/9vHxGTx4\nsHY1ICAgKirqX//61/Dhwz0w8bdcmTw4OPjuu+/WLtfV1W3dujUqKio6OtoD436XK8O7sqYt\ndo56n4aGhpCQkJaWlkcffXT+/PnNzc0LFiyor6+/6OLW1ta33nprypQpfeFTWKdPnxaRrKys\n5OTk119/PSEh4Q9/+ENNTU3bNbW1tYGBgW0/bhMcHNxuTe9zZXJNa2vr5MmT77vvvuLi4tTU\nVB8fn14ftj1Xhnf9CWoIh/cJDg7esGHD7Nmzhw4dOnTo0OTkZKvV+tlnn1108e7du61W6/jx\n43t5yE4kJiZGRUUFBgY++OCDOp1u37597Rb02Q/pdTm5iOj1+j/96U/PP/98UFBQSkrKuXPn\nen/Oi3JleFfWaAiH1/P39+/fv/+pU6cu+tXs7OyxY8f2kc9uhoWFybcfkRARg8EQFhZWXV3d\ndk1ISEhtba2jzQnNNTU1Hv8zOa5M7hQVFfWjH/0oOTm5pqbm008/7b0pO+DK8EpPUAiHNyoq\nKlq1alVLS4t21Wq1VlVVDRw48MKV9fX1ubm5o0eP7t0BOxQWFhYaGnr48GHtqnboZ8CAAW3X\nDBkyxGazHTt2TLtaW1tbXFx81VVX9fas3+XK5Lm5uQ899FBTU5N2VafT9YXNQ3FteFfWtNUn\nnhiUhIWFff755y0tLVOmTLHb7Rs2bAgICBg7dqyIfPTRR1ardeLEidrKo0eP2u32K664wqPz\n/pter584ceLGjRujoqKioqIyMjJMJpN2soBz8rCwsDFjxrzyyitPPPGEr6/vX//61yuvvPLq\nq6/u+5MPGTLEarX+8Y9/nDp1qo+Pz3vvvWe1Wq+77jrPTu7i8J2suSg+5OaVjh8/vm7dOu0A\nREJCwsyZM7VfDsuWLautrU1NTdWWffLJJ2lpaZs3b+4jv/pEpLW1NT09PSsr69y5cwkJCY8+\n+qh23KHt5A0NDWvXrs3NzbXb7T/84Q8ffvhhj2+qiGuTFxUVrVu37uDBgzqdLiYm5te//rVn\nDwY5uTJ8R2suinAAUMY+DgDKCAcAZYQDgDLCAUAZ4QCgjHAAUEY4ACgjHACUEQ4AyggHAGWE\nA4Cy/wNAWghinHn++wAAAABJRU5ErkJggg=="
          },
          "metadata": {
            "image/png": {
              "width": 180,
              "height": 180
            }
          }
        }
      ]
    }
  ]
}