Module5-BayesianThinkingAndMCMC/M5A_Bayesian_Thinking.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Bayesian Thinking\n",
    "<html>\n",
    "    <summary></summary>\n",
    "         <div> <p></p> </div>\n",
    "         <div style=\"font-size: 20px; width: 800px;\"> \n",
    "              <h1>\n",
    "               <left>Bayesian Thinking</left>\n",
    "              </h1>\n",
    "              <p><left>============================================================================</left> </p>\n",
    "<pre>Course: BIOM 421, Spring 2024\n",
    "Instructor: Brian Munsky\n",
    "Authors: Drs. Kaan Öcal, Huy Vo, Brian Munsky\n",
    "Contact Info: munsky@colostate.edu\n",
    "</pre>\n",
    "         </div>\n",
    "    </p>\n",
    "\n",
    "</html>\n",
    "\n",
    "<details>\n",
    "  <summary>Copyright info</summary>\n",
    "\n",
    "```\n",
    "Copyright 2024 Brian Munsky\n",
    "\n",
    "Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:\n",
    "\n",
    "1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.\n",
    "\n",
    "2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.\n",
    "\n",
    "3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.\n",
    "\n",
    "THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS \"AS IS\" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n",
    "```\n",
    "<details>\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "kfORGQnMZfb4"
   },
   "outputs": [],
   "source": [
    "# Importing packages\n",
    "import numpy as np\n",
    "import scipy as sp\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "kGUws9EGuUuI"
   },
   "source": [
    "## Introduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Vl_sATEldNw4"
   },
   "source": [
    "Bayesian statistics is a cornerstone of modern statistics. In this lecture we will learn what it is, and where it can be used. Before tackling this notebook, we highly recommend the video [\"Bayes theorem, the geometry of changing beliefs\"](https://www.youtube.com/watch?v=HZGCoVF3YvM) by 3Blue1Brown."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "lEP3hqGVuaQy"
   },
   "source": [
    "## Why Bayesian inference?\n",
    "\n",
    "Biological experiments are often slow and costly, and due to the delicate nature of living cells we cannot directly measure many of the things we care about. When we do manage to measure them, our measurements are typically subject to a lot of noise - we saw this in our analysis of smFISH microscopy images. If we want to study biology, therefore, we need to be careful about the conclusions we draw given imperfect data, and Bayesian statistics provides a general framework that allows us to deal with the noise and uncertainty we encounter in practice.\n",
    "\n",
    "The workhorse of Bayesian inference is **Bayes' Theorem**:\n",
    "\n",
    "$$ p(x \\, | \\, \\mathcal{Data}) = \\frac {p(\\mathcal{Data} \\, | \\, x) \\, p(x)}{p(\\mathcal{Data})} $$\n",
    "\n",
    "Here $x$ is quantity that we wish to infer, such as the bias of a coin, or the population of deer in Lory State Park. $\\mathcal{Data}$ is what's observed, which could be the outcome of repeated coin flips, or reported deer sightings in the past year. Bayes' Theorem tells us how observing $\\mathcal{Data}$ affects our knowledge of the unknown quantity $x$, based on our prior knowledge of $x$ (represented by the prior distribution $p(x)$), and the likelihood function $p(\\mathcal{Data} \\, | \\, x)$. \n",
    "\n",
    "The normalising constant in the denominator of Bayes' Theorem is required to ensure that the posterior probabilities over $x$ sum to $1$. The posterior is meaningless without this. The normalising constant can be explicitly written as\n",
    "\n",
    "$$ p(\\mathcal{Data}) = \\int p(\\mathcal{Data} \\, | \\, x) \\, p(x) \\, dx $$\n",
    "\n",
    "The normalization constant is the total probability of observing $\\mathcal{Data}$ summed over *all* possible values of $x$. In order for us to know whether a specific value of $x$ is a good explanation of our data, we need to consider how well that explanation performs *relative to all other explanations*.\n",
    "\n",
    "To summarize the elements of the Bayes' Theorem (remember and understand this!!)\n",
    "* $x$ = quantity that we wish to infer,\n",
    "* $p(x)$ = our **prior** knowledge of $x$ before doing the experiment.\n",
    "* $\\mathcal{Data}$ = our data or **evidence**.\n",
    "* $p(\\mathcal{Data} \\, | \\, x)$ = our model predictions or **likelihood function**.\n",
    "* $p(\\mathcal{Data})$ = the **normalization constant**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example 1 - Librarians and Farmers\n",
    "In [3Blue1Brown's video on Bayes' theorem](https://www.youtube.com/watch?v=HZGCoVF3YvM) you will meet Steve, and you will be given four facts:\n",
    " \n",
    " 1) Steve is shy and withdrawn, \n",
    " 2) 40% of librarians are shy and withdrawn\n",
    " 3) 10% of farmers are shy and withdrawn\n",
    " 4) There are 20 times as many farmers as there are librarians\n",
    " 5) (assumed) Steve is known to be either a farmer or a librarian.\n",
    "\n",
    "* How would you guess if Steve is a librarian or a farmer?\n",
    "* What do you think the probability is that Steve is a librarian?\n",
    "\n",
    "How do you rationally answer these questions?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's define the variables of our problem:\n",
    "P_data_given_librarian = 0.4\n",
    "P_data_given_farmer = 0.6\n",
    "P_prior_librarian = 1 / (1 + 20)\n",
    "P_prior_farmer = 1 - P_prior_librarian\n",
    "\n",
    "# Now we can calculate the probability of the data given the model:\n",
    "P_shy_and_withdrawn = P_data_given_librarian*P_prior_librarian + P_data_given_farmer*P_prior_farmer\n",
    "\n",
    "# And the probability of the model given the data:\n",
    "P_librarian_given_data = P_data_given_librarian*P_prior_librarian / P_shy_and_withdrawn\n",
    "P_farmer_given_data = P_data_given_farmer*P_prior_farmer / P_shy_and_withdrawn\n",
    "\n",
    "print(f\"The probability that Steve is a librarian is {P_librarian_given_data:.2f}\")\n",
    "print(f\"The probability that Steve is a farmer is {P_farmer_given_data:.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example 2 - COVID Test\n",
    "According to [Cochrane](https://www.cochrane.org/CD013705/INFECTN_how-accurate-are-rapid-antigen-tests-diagnosing-covid-19) I found that at home rapid antigen tests have the following statistics:\n",
    " \n",
    " 1) The probability of a false negative is 18%.\n",
    " 2) The probability of a true positive is 92%.\n",
    " 3) The probability of a false positive is 0.4%.\n",
    " 4) The probability of a true negative is 99.6%.\n",
    "\n",
    "We want to know:\n",
    "\n",
    "* If your test is positive, what is the probability that you have COVID?\n",
    "* If your test is negative, what is the probability that you have COVID?\n",
    "\n",
    "It might not be clear at first, but to answer this question you need to know that about 5% of people are infected at the current time."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "|| COVID | No Covid |\n",
    "| -: | --- | -: |\n",
    "|Positive Test | 92% | 0.4% |\n",
    "|Negative Test | 18% | 99.6% | \n",
    "||||\n",
    "|Total|100%|100%|"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the variables of our problem:\n",
    "P_positive_given_covid = 0.92\n",
    "P_negative_given_covid = 0.08\n",
    "P_positive_given_no_covid = 0.08\n",
    "P_negative_given_no_covid = 0.92\n",
    "P_prior_covid = 0.05\n",
    "P_prior_no_covid = 1 - P_prior_covid\n",
    "\n",
    "# Calculate the probability of the data given the model:\n",
    "P_positive = P_positive_given_covid*P_prior_covid + P_positive_given_no_covid*P_prior_no_covid\n",
    "P_negative = P_negative_given_covid*P_prior_covid + P_negative_given_no_covid*P_prior_no_covid\n",
    "\n",
    "# Calculate the probability of the model given the data:\n",
    "P_covid_given_positive = P_positive_given_covid*P_prior_covid/P_positive\n",
    "P_covid_given_negative = P_negative_given_covid*P_prior_covid/P_negative\n",
    "\n",
    "print(f\"The probability that a person has Covid given a positive test is {P_covid_given_positive:.3f}\")\n",
    "print(f\"The probability that a person has Covid given a negative test is {P_covid_given_negative:.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "vxylSMJCnjcN"
   },
   "source": [
    "## **Example 3:** Flipping coins\n",
    "\n",
    "*In this example, we are going to move away from the binary (yes/no) problem and try to estimate a parameter on the real line.*\n",
    "\n",
    "We can visualise the Bayesian approach to parameter estimation using a simple coin example. Suppose I have a coin which has a fixed probability $h$ to land on heads each time it is flipped. How much can we say about $h$ by observing how the coin behaves? We flip the coin $8$ times and observe the sequence $H, H, H, H, T, H, H, H$.\n",
    "\n",
    "Let us analyse the individual parts in Bayes' Theorem for this example:\n",
    "1. **The prior $p(h)$:** We assume the coin is fair, but allow for the possibility of a biased coin. Our prior should be peaked at $h = \\frac 1 2$, but allow all values between $0$ and $1$. A convenient probability distribution for numbers between $0$ and $1$ is the [Beta distribution](https://en.wikipedia.org/wiki/Beta_distribution):\n",
    "\n",
    "$$ p(h) = \\textrm{Beta}(h; 4, 4) $$\n",
    "\n",
    "2. **The likelihood $p(\\mathcal{Data}  \\, | \\, h)$:** Each flip has a probability $h$ to land on heads, independently of the others. We can write down the likelihood explicitly:\n",
    "$$ p(\\mathcal{Data} \\, | \\, h) = h^7 (1-h) $$\n",
    "Please note: this is the likelihood of the sequence in its known order. If the data were unordered (i.e., 7 heads and one tails but in an unknown order), you would need to use the Binomial distribution. Since the the order of the coin flips does not matter here, both approaches give the same posterior.\n",
    "\n",
    "3. **The posterior $p(h \\, | \\, \\mathcal{Data})$:** In this example, a little bit of maths allows us to compute this explicitly. We will later see how to do this using MCMC.\n",
    "\n",
    "In the exercises you will be asked to compute the posterior, which requires knowing the normalisation constant $p(\\mathcal{Data})$. Even for this simple example, computing the normalisation constant is not trivial!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "zcoQRBD6ivh5"
   },
   "outputs": [],
   "source": [
    "# Estimating the fairness of a coin based on repeated coin flips\n",
    "data_ct = [ 1, 1, 1, 1, 0, 1, 1 ]  # 1 for heads, 0 for tails\n",
    "n_heads = np.sum(data_ct)\n",
    "\n",
    "# Prior for the coin toss example\n",
    "def prior_ct(h):\n",
    "  return sp.stats.beta(4,4).pdf(h)\n",
    "\n",
    "# Likelihood function\n",
    "def likelihood_ct(data, h):\n",
    "  n_heads = sum(data)\n",
    "  n_tails = len(data) - n_heads\n",
    "\n",
    "  return h ** n_heads * (1-h) ** n_tails\n",
    "\n",
    "# Posterior\n",
    "# We will ignore the normalisation constant for now; it does not appear in the plot\n",
    "def posterior_unnormalised_ct(h, data):\n",
    "  return likelihood_ct(data, h) * prior_ct(h)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 360
    },
    "id": "h68a5Yk9xWOB",
    "outputId": "cf7c595c-ecbd-4052-86eb-8b173c7a8cb1"
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, sharex=True, figsize=(10,3))\n",
    "\n",
    "xx = np.arange(0, 1, 0.002)\n",
    "\n",
    "# Plot prior\n",
    "axes[0].plot(xx, prior_ct(xx), color=\"blue\")\n",
    "axes[0].fill_between(xx, 0, prior_ct(xx), color=\"blue\", alpha=0.3)\n",
    "axes[0].set_title(\"Prior\\n\")\n",
    "axes[0].set_xlabel(\"Head probability $h$\")\n",
    "axes[0].set_yticks([])\n",
    "axes[0].set_ylim(bottom=0)\n",
    "\n",
    "# Plot likelihood\n",
    "axes[1].plot(xx, likelihood_ct(data_ct, xx), color=\"orange\")\n",
    "axes[1].set_title(\"Likelihood\\n\")\n",
    "axes[1].set_xlabel(\"Head probability $h$\")\n",
    "axes[1].set_yticks([])\n",
    "axes[1].set_ylim(bottom=0)\n",
    "\n",
    "# Plot the (unnormalised) posterior\n",
    "axes[2].plot(xx, posterior_unnormalised_ct(xx, data_ct), color=\"green\")\n",
    "axes[2].fill_between(xx, 0, posterior_unnormalised_ct(xx, data_ct), color=\"green\", alpha=0.3)\n",
    "axes[2].set_title(\"Posterior\\n($\\\\propto$Prior $\\\\times$ Likelihood)\")\n",
    "axes[2].set_xlabel(\"Head probability $h$\")\n",
    "axes[2].set_yticks([])\n",
    "axes[2].set_ylim(bottom=0)\n",
    "axes[2].set_xlim(0, 1);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "kX2_CDnUISQm"
   },
   "source": [
    "## **Example 4:** Gene expression\n",
    "\n",
    "Consider a bacterial gene $G$ that produces mRNA $M$ at a fixed but unknown rate $\\sigma$. mRNA typically gets degraded at a fixed rate $d$ that can be measured directly; in this example, the degradation rate is measured to be about $d = 0.5/\\textrm{h}$.\n",
    "\n",
    "$$ G \\longrightarrow G + M \\qquad\\qquad M \\rightarrow \\emptyset $$\n",
    "\n",
    "We want to infer the production rate $\\sigma$ by measuring gene expression in 30 cells using an smFISH experiment. After the experiments have been performed, the data looks like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "6xNcFz_oxvgd"
   },
   "outputs": [],
   "source": [
    "data_ge = [ 20, 20, 19, 16, 15, 22, 17, 27, 17, 17, 21, 21, 16, 22, 25, 22, 23,\n",
    "            20, 21, 16, 16, 18, 16, 21, 17, 21, 25, 16, 15, 23 ]\n",
    "\n",
    "deg = 0.5     # degradation rate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 311
    },
    "id": "t06tiXQ4Juaq",
    "jupyter": {
     "source_hidden": true
    },
    "outputId": "56f74faf-2949-4565-b01a-59314c26cfba"
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(3,3))\n",
    "\n",
    "ax.hist(data_ge, bins=21, range=(9.5, 30.5))\n",
    "\n",
    "ax.set_xlabel(\"mRNA\")\n",
    "ax.set_ylabel(\"Frequency\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "wtPva-04K500"
   },
   "source": [
    "Let us review the individual parts in Bayes' Theorem for this example.\n",
    "1. **The prior $p(\\sigma)$:** Bacteria tend to be small and express only few mRNA molecules. A realistic production rate could be about $5/\\mathrm{h}$, but this depends on the gene. For our prior we choose an exponential distribution with mean $5$:\n",
    "\n",
    "$$ p(\\sigma) = \\textrm{Exp}(\\sigma; 1/5) $$\n",
    "\n",
    "2. **The likelihood $p(\\mathcal{Data}  \\, | \\, \\sigma)$:** The cells have been kept in the lab for a while, so we assume that mRNA production has reached steady state. At steady state, mRNA counts for this model are Poisson distributed with mean $\\sigma / d$:\n",
    "$$ p(n \\, | \\, \\sigma) = \\textrm{Poisson}(n; \\sigma / d) $$\n",
    "  Here $n$ is the number of mRNA measured. The total likelihood is just the product of the likelihoods for each individual cell, since we assume that each cell is independent of the others:\n",
    "\n",
    "$$ p(\\mathcal{Data} \\, | \\, \\sigma) = \\prod_{i=1}^{30} \\textrm{Poisson}(n_i; \\sigma / d) $$\n",
    "\n",
    "3. **The posterior $p(\\sigma \\, | \\, \\mathcal{Data})$:** In this example we can still compute this explicitly. This won't be possible for more complex models of gene expression. The difficulty here is again the normalisation constant.\n",
    "\n",
    "**Important:** In the Bayesian setting we have to decide the prior *before* looking at the data - the data only comes in via Bayes' formula."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "b-ZFr8NoyYvv"
   },
   "outputs": [],
   "source": [
    "# Estimating the mRNA production rate based on snapshot data\n",
    "# Prior for the gene expression model\n",
    "def prior_ge(h):\n",
    "  return sp.stats.gamma(1, scale=5).pdf(h)\n",
    "\n",
    "# Likelihood\n",
    "def likelihood_ge(data, sigma):\n",
    "  ret = 1\n",
    "  for n in data:\n",
    "    ret *= sp.stats.poisson(sigma / deg).pmf(n)\n",
    "\n",
    "  return ret\n",
    "\n",
    "# Posterior (unnormalised, again)\n",
    "def posterior_unnormalised_ge(h, data):\n",
    "  return likelihood_ge(data, h) * prior_ge(h)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 357
    },
    "id": "HnV559bxJm2h",
    "jupyter": {
     "source_hidden": true
    },
    "outputId": "47060302-e442-4250-cee7-784c0f16983d"
   },
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 3, sharex=True, figsize=(10,3))\n",
    "\n",
    "xx = np.arange(0, 20, 0.01)\n",
    "\n",
    "# Plot the prior\n",
    "axes[0].plot(xx, prior_ge(xx), color=\"blue\")\n",
    "axes[0].fill_between(xx, 0, prior_ge(xx), color=\"blue\", alpha=0.3)\n",
    "axes[0].set_title(\"Prior\\n\")\n",
    "axes[0].set_xlabel(\"Production rate $\\\\sigma$\")\n",
    "axes[0].set_yticks([])\n",
    "axes[0].set_ylim(bottom=0)\n",
    "\n",
    "# Plot the likelihood\n",
    "axes[1].plot(xx, likelihood_ge(data_ge, xx), color=\"orange\")\n",
    "axes[1].set_title(\"Likelihood\\n\")\n",
    "axes[1].set_xlabel(\"Production rate $\\\\sigma$\")\n",
    "axes[1].set_yticks([])\n",
    "axes[1].set_ylim(bottom=0)\n",
    "\n",
    "# Plot the posterior\n",
    "axes[2].plot(xx, posterior_unnormalised_ge(xx, data_ge), color=\"green\")\n",
    "axes[2].fill_between(xx, 0, posterior_unnormalised_ge(xx, data_ge), color=\"green\", alpha=0.3)\n",
    "axes[2].set_title(\"Posterior\\n($\\\\propto$ Prior $\\\\times$ Likelihood)\")\n",
    "axes[2].set_xlabel(\"Production rate $\\\\sigma$\")\n",
    "axes[2].set_yticks([])\n",
    "axes[2].set_ylim(bottom=0)\n",
    "axes[2].set_xlim(0, 20);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "_ncywF5MQqVw"
   },
   "source": [
    "We see that the posterior is tightly peaked around $10$, which suggests that we have enough data to make confident estimates of $\\sigma \\approx 10 \\pm 1$.\n",
    "\n",
    "**Note:** Bayesian inference does not produce confidence intervals. Credible intervals take the place of confidence intervals in Bayesian analysis."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "CYzb62japVFG"
   },
   "source": [
    "<!-- This simple example showcases the main aspects of Bayesian inference: choosing a prior, computing the likelihood, and computing the posterior.  -->\n",
    "\n",
    "## Bayesian inference in practice\n",
    "\n",
    "In both examples we could use Bayes' theorem to explicitly compute the posteriors (see exercises below). In practice we often encounter more complex models, and each of the three terms in Bayes' formula can cause headaches.\n",
    "1. **The prior $p(x)$:** We usually don't know much about $x$ before collecting data, but we need to choose a prior distribution. Sometimes it's not obvious what prior is \"correct\", and choosing the wrong prior can significantly bias our results if we don't have enough data.\n",
    "2. **The likelihood $p(\\mathcal{Data} \\, | \\, x)$:** This is often difficult to compute, depending on the model and the type of data. Some models have such complicated likelihoods that we cannot feasibly compute them at all! In the next lecture we will encounter likelihoods that are based on the Chemical Master Equation, which is solvable depending on the system.\n",
    "3. **The posterior $p(x \\, | \\, \\mathcal{Data})$:** This requires the normalising constant $p(\\mathcal{Data})$, which in turn requires computing the likelihood for *all* $x$. This is often the hardest part of Bayesian inference."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Jq6Ashz3t14p"
   },
   "source": [
    "## Exercises\n",
    "<!-- 1. Using the formula above, compute the mean and variance of the Gaussian distribution. Can you see what the factor of $1/\\sqrt{2\\pi}$ is doing in the formula?\n",
    "\n",
    "  (*Hint:* For the mean, you can bypass explicit computations by noticing a certain symmetry. For the variance, use integration by parts.) -->\n",
    "\n",
    "1. Review the explicit formula we gave for the normalising constant:\n",
    "$$ p(\\mathcal{Data} \\, | \\, x) = \\int p(\\mathcal{Data} \\, | \\, x) \\, p(x) \\, dx $$\n",
    "Verify that this choise ensures that the posterior is in fact a probability distribution, that is,\n",
    "\n",
    "$$ \\int p(x \\, | \\, \\mathcal{Data}) \\, dx = 1 $$\n",
    "\n",
    "2. **a.** The prior distribution for the coin flip example has the following formula:\n",
    "\n",
    "  $$ \\textrm{Beta}(h; 4, 4) = \\frac{7!}{3! \\times 3!} h^3 (1-h)^3 $$\n",
    "\n",
    "  for $h$ between $0$ and $1$. Can you verify that this distribution is normalised correctly, ie. that it integrates to $1$?\n",
    "\n",
    "  (*Hint:* Use integration by parts. As often as necessary.)\n",
    "\n",
    "  **b.** (harder) Can you compute the normalising constant $p(\\mathcal{Data})$ explicitly?\n",
    "\n",
    "  (*Hint:* This is similar to **a**, but involves even more integration by parts.)\n",
    "\n",
    "3. Can you compute the normalising constant $p(\\mathcal{Data})$ for the gene expression example?\n",
    "\n",
    "  (*Hint:* You can use the integral identity\n",
    "\n",
    "  $$ \\int_0^\\infty x^m \\exp(-a x) \\, dx = \\frac{m!}{a^{m+1}}, $$\n",
    "\n",
    "  which of course can be verified using integration by parts.)"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}