Ctrl+K

Lab 9: Bayesian Analysis

Lab 9: Bayesian Analysis#

For the class on Monday, February 12th

Tip

Looking for hints?

A. Comparing Samples#

We have two sets of data points. Data Set 1 has \(N_X\) points \((X_1, X_2, ..., X_{N_X})\). Data Set 2 has \(N_Y\) points \((Y_1, Y_2, ..., Y_{N_Y})\).

We know both are sampled from an exponential distribution, but potentially with different values of \(\lambda\). We want to test if these two data sets come from the same exponential distribution (i.e., same \(\lambda\)).

Let’s say Data Set 1 is sampled from \(\text{Exp}(\lambda_X=a+b)\) and Data Set 2 from \(\text{Exp}(\lambda_Y=a-b)\).

  1. Write down the likelihood function in terms of \(a\), \(b\), \({X_i}\) and \({Y_i}\).

  2. Let’s take \(\text{Gamma}(\alpha=2, \beta=2) = 4ae^{-2a}\) as the prior for \(a\), and \(\mathcal{N}(\mu=0, \sigma=1) = \frac{1}{\sqrt{2\pi}}e^{-b^2/2}\) as the prior for \(b\). Write down the posterior function (up to a normalizing constant) in terms of \(a\), \(b\), \({X_i}\) and \({Y_i}\).

// Write your answers to Part A here

Demo: Find confidence interval#

The code below demonstrates how to numerically find the 95% confidence interval around 0 for the standard normal distribution.

import numpy as np
from scipy.integrate import cumulative_trapezoid
from scipy.stats import norm

x = np.linspace(0, 3, 30001)
y = np.exp(-0.5*x*x) / np.sqrt(2*np.pi)
y_int = cumulative_trapezoid(y, x=x, initial=0)
print("calculated value:", x[np.searchsorted(y_int, 0.95/2)])
print("      true value:", norm.isf((1-0.95)/2))

calculated value: 1.9600000000000002
      true value: 1.959963984540054

B. Coin Flips#

Assume this coin shows heads with a probability of \(p\), and tails with a probability of \(1-p\).

B1. A simple case#

We flipped the coin \(N=100\) times and got \(K=60\) heads.

  1. Plot the posterior of \(p\) assuming a uniform prior.

  2. What is the maximum a posteriori probability (MAP) estimate?

  3. What is the 95% confidence interval of your estimate of \(p\)?

  4. Is your result from (3) consistent with what you obtained in Lab 3, Part A?

// Write your answers to Part B1 here

# Include your implementation for Part B1 here

B2. Small number of flips#

We flipped a coin \(N=3\) times and got \(K=0\) heads (no heads).

  1. Plot the posterior of \(p\) assuming a uniform prior.

  2. What is the maximum a posteriori probability (MAP) estimate?

  3. What is the 95% confidence interval of your estimate of \(p\)?

  4. Make a prediction for the number of heads you would obtain if you flip this particular coin 100 times.

// Write your answers to Part B2 here

# Include your implementation for Part B2 here

B3. Impact of the prior#

  1. Repeat B1 and B2 (all steps in both sections), but replace the prior with \(\text{Beta}(\alpha=2, \beta=2) = 6 p(1-p)\).

  2. Do your answers with this new prior agree with what you expect? Briefly explain.

// Write your answers to Part B3 here

# Include your implementation for Part B3 here

Final question#

Roughly how much time did you spend on this lab (not including the time you spent in class)?

// Write your answers here

Tip

How to submit this notebook on Canvas?

Once you complete this lab and have all your code, text answers, and all the desired results after running the notebook properly displayed in this notebook, please convert the notebook into HTML format by running the following:

jupyter nbconvert --to html /path/to/labs/09.ipynb

Then, upload the resulting HTML file to Canvas for the corresponding assignment.

Contents