<jupyter_start><jupyter_text>***Your Name Here***<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
%matplotlib inline
from sklearn.neighbors import KernelDensity<jupyter_output><empty_output><jupyter_text>Prelab 9 (Last!)***This exercise builds on the ideas that you encountered in Lab 8, so you should complete Lab 8 before beginning it.*** Let's be explicit about what you did in Lab 8, namely:   1. You ***measured*** the means of two arrays of values (Hubble's distances and recessional velocities) and computed their covariance matrix (which captures both the standard deviations of each quantity and the correlation coefficient of their relationship to one another).  2. You used these computed values to create a statistical ***model*** of the data (a bivariate gaussian with parameters $\mu_d, \mu_v, \sigma_d, \sigma_v, r$ derived from the data).   3. You ***synthesized samples*** from this data model in a Monte Carlo simulation and computed two statistics with real physical meaning, namely: (a) the correlation coefficient r, which reveals the strength of the relationship between recessional velocity and distance, and (b) the slope of best fit line (also known as the Hubble constant $H_0$) which describes how a change in one variable affects the other, on average.  4. You explored the variation in this ***synthesized sampling distribution*** by plotting histograms of the best fit values for r and $H_0$ from all of your random draws of 24 galaxies. In a Bayesian framework, we might say that these histograms are visualizations of the Likelihood, or a way of looking at the strength of the evidence in favor of the experimental results (Hubble's data) under the (bivariate gaussian) model. We could turn this into a posterior probability by multiplying it by some prior and dividing by a normalization constant (because probabilities always sum to 1). Indeed, this is what you will do below. 5. You then completed a classical hypothesis test. At the end of this prelab, you will reimagine that as a Bayesian hypothesis test (should you choose to do so in Part E).  Part AFirst, bring over the code that you used to complete steps 1-4 of Lab 8. Add a comment to the top of each cell describing the big picture of what it contains.  Part BBelow these instructions is a cell that makes a "kernel density estimate", which is a way of turning your histogram into a continuous function with a defined value at a range of x values. You can use this smooth function as a Likelihood function in a Bayesian framework. Comment each line in the provided code and play around with the tunable parameters (especially the bandwidth and kernel keywords). Once you feel you understand what the code is doing in general, describe it in the cell provided, including the effects of changing the values of the tunable parameters. Use visualizations to support your arguments.   Part C1. Compute the area under the entire KDE by summing the appropriate array. What is it equal to and why?   2. Write a function that will sum under the KDE for an arbitrary range of x values and return that sum. Note the input should be x values, ***not*** indices, so you will need to translate values (e.g. $x_1$=100, $x_2$=300) to their indices (or nearest index, e.g. i=0, j=200) within the discrete x array so that you can sum over the correct range of indices in the prob array. Note that what you're doing here is essentially computing the value of the finite integral $$prob = norm\int_{x_1}^{x_2}L(x)dx$$ where L(x) is the likelihood, prob is the probabilty of obtaining a value between $x_1$ and $x_2$, and norm is a normalization factor, where necessary. Part DNow let's add in some priors.1. First, devise an ***uninformative prior*** over the same range of x values as the KDE and plot it.    *Hint: some manipulation may be required to make sure your prior is a plot of prior probabilities. What special property do probabilties have under summing?)*  2. Now, let's develop two different ***informative priors*** according to the following scenarios:    a) A similar sample of galaxy distances and recessional velocities compiled by Dr. <NAME> is uncovered. Let's not assume anything crazy about these data - just that he also measured 24 random galaxies drawn from roughly the same bivariate gaussian distribution as we assumed Hubble's came from (so basically let's assume the Hubble model is a correct representation of the true universe). Of course, in his case, he has drawn a different random realization of 24 points from this distribution. To simulate this, pick one set of distances and recessional velocities from your Monte Carlo simulator that has an $H_0$ value a little different from Hubble's. Use this information to create a prior for $H_0$ assuming Dr. Bubble's data is all that we know about. Use the values of recessional velocity and distance for your chosen random realization as a starting point, and the code from Part A to draw samples from a slightly different bivariate gaussian centered on the Bubble data values. Use it to make a histogram and turn this into a Likelihood, which you will use as a prior.        b) Another mysery scientist (you can name her if you like) has completed a comparable study that is MORE informative than Hubble's data. It has roughly the same average values for recessional velocities and distances, but the sampling distribution is narrower. (how can you accomplish this? There are several ways!). In this case, let's assume our Hubble data are the prior and these data are the Likelihood. Describe the technique you used to make this prior more informative than the Hubble data and how you know you acheived your goal.   *Tips:*    * *In most cases, you will want to use the specifications for the informative prior to create a new sampling distribution histogram using your code in A, then convert it to a kde then likelihood as you did in part B*   * *If creating a new covariance matrix, remind yourself of what the values in the covariance matrix represent and remember that the matrix is symmetric on either side of the diagonal, so you may have to change two values instead of just one)*    * *You may wish to show your histograms for the original (Hubble) and new Monte Carlo models side by side to be sure that you understand the results of your modifications*  3. For each of these priors (1, 2a, and 2b), make a plot showing the prior, likelihood, and posterior (derived from the product of the prior and the likelihood, but with a little additional manipulation - see if you can figure out what is missing on your own. Bayes' theorem may provide some hints) on the same plot. The plot should include a legend and axis labels. You may wish to use fill_beteween as in the example below and an alpha (transparency) value so that you can visualize overlap between the three curves via combined colors (e.g. the overlap region between a red likelihood and a blue posterior is a purple color)  4. In a markdown cell, summarize the "takeaway messages" of the three plots you made in Part 3. What do they reveal about priors, likelihoods, and posteriors?   [Optional - Challenge Exercise] Part ECalculate Bayes' factor for the following two hypotheses, implementing one of the informative priors from the example above:   1. The true value of $H_0$ is within 1 standard deviation of the measured value.  2. The true value of $H_0$ is NOT within 1 standard deviation of the measured value.<jupyter_code>## Part A code placeholder. Make sure to comment each cell with the "big picture"
#Part B - sample code to comment
#hlist is the array of $H_0$ values returned by your Monte Carlo simulation
x_d=np.arange(200,700)
kde = KernelDensity(bandwidth=10.0, kernel='gaussian')
kde.fit(hlist[:, None])
logprob = kde.score_samples(x_d[:, None])
prob=np.exp(logprob)
plt.fill_between(x_d, prob, alpha=0.5)
#Part B exploration (copy and paste the cell above and manipulate it to explore what it's doing)<jupyter_output><empty_output><jupyter_text>Part B explanation. Make sure to use plots to support your arguments.<jupyter_code>#Part C #1 Code<jupyter_output><empty_output><jupyter_text>Part C 1 Explanation<jupyter_code>#Part C #2 Function definition - make sure to comment
#Part C #2 Function test statements
#Part D #1 - Define the uninformative prior
#Part D #1 - Plot the uninformative prior
#Part D #2 - Define the first informative prior (Edgar Bubble data)
#Part D #2 - Define the second informative prior (more informative)<jupyter_output><empty_output><jupyter_text>Part D 2 - informative prior 2 explanation<jupyter_code>#Part D #2 (optional but highly suggested) - visualize Hubble, Bubble, and more informative sampling histograms side by side
#Part D #3 - Plot of prior, likelihood, and posterior probability for uninformative prior
#Part D #3 - Plot of prior, likelihood, and posterior probability for informative prior #1
#Part D #3 - Plot of prior, likelihood, and posterior probability for informative prior #2<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup