{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# SMIPP 21/21 - Exercise Sheet 3\n", "\n", "## Prof. Dr. K. Reygers, Dr. R. Stamen, Dr. M. Völkl\n", "\n", "## Hand in by: Thursday, November 11th: 12:00\n", "### Submit the file(s) through the Übungsgruppenverwaltung\n", "\n", "\n", "### Names (up to two):\n", "### Points: " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3.1 Error propagation for correlated variables (10 points)\n", "\n", "a) Propagate uncertainties for the following expressions using [SymPy](https://www.sympy.org) following the [example notebook](https://nbviewer.jupyter.org/urls/www.physi.uni-heidelberg.de/~reygers/lectures/2020/smipp/gaussian_error_propagation_correlated_variables.ipynb) of the lecture. Assume that the variables are correlated and that the correlation coefficient is $\\rho$.\n", "\n", "i) The acceleration of gravity with a simple pendulum is given by the following formula:\n", "$$g = 4 \\pi^2 \\frac{L}{T^2}$$\n", "The relevant variables are the length $L$ of the pendulum and the period $T$ with the corresponding errors $\\sigma_L$ and $\\sigma_T$.\n", "\n", "ii) The determination of various important physics quantities in particle physics is based on the measurement of asymmetries. One of these asymmetries is the so called forward-backward asymmetry which is defined as\n", "$$A = \\frac{F-B}{F+B}$$\n", "$F$ is the number of certain particles which move into a predefined forward direction and $B$ is the number of particles moving into the corresponding backward direcion with the errors $\\sigma_F$ and $\\sigma_B$.\n", "\n", "b) The scattering angle and the radial distance of a certain particle can be determined from a position measurement $(x,y)$ \n", "$$r = \\sqrt{x^2 + y^2}, \\quad \\theta = \\mathrm{atan2}(y, x)$$\n", "You find more on the [atan2](https://en.wikipedia.org/wiki/Atan2) function on wikipedia. The position ($x$,$y$) is measured with the corresponding uncertainties $\\sigma_x$ and $\\sigma_y$. Write a python function that returns the covariance matrix $U$ of $r$ and $\\theta$ for a given covariance matrix $V$ of $x$ and $y$. Determine $U$ under the assumption that $x$ and $y$ are uncorrelated. Hint: use the formula on page 17 of the slides of chapter 3.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solution" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Your Solution" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3.2 Covariance and correlation coefficient (10 points)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "a) Plot the four data sets ($y$ vs. $x$)\n", "\n", "b) Calculate the means $\\mu_x$, $\\mu_y$, the unbiased sample variances $\\sigma_x^2$, $\\sigma_y^2$, the unbiased covariance, and the corresponding correlation coefficient for the four data sets below without using `numpy.mean`, `numpy.var`, `numpy.cov` and `numpy.corrcoef`. Are you surprised?" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "x1 = np.array([10.0, 8.0, 13.0, 9.0, 11.0, 14.0, 6.0, 4.0, 12.0, 7.0, 5.0])\n", "y1 = np.array([8.04, 6.95, 7.58, 8.81, 8.33, 9.96, 7.24, 4.26, 10.84, 4.82, 5.68])\n", "\n", "x2 = np.array([10.0, 8.0, 13.0, 9.0, 11.0, 14.0, 6.0, 4.0, 12.0, 7.0, 5.0])\n", "y2 = np.array([9.14, 8.14, 8.74, 8.77, 9.26, 8.10, 6.13, 3.10, 9.13, 7.26, 4.74])\n", "\n", "x3 = np.array([10.0, 8.0, 13.0, 9.0, 11.0, 14.0, 6.0, 4.0, 12.0, 7.0, 5.0])\n", "y3 = np.array([7.46, 6.77, 12.74, 7.11, 7.81, 8.84, 6.08, 5.39, 8.15, 6.42, 5.73])\n", "\n", "x4 = np.array([8.0, 8.0, 8.0, 8.0, 8.0, 8.0, 8.0, 19.0, 8.0, 8.0, 8.0])\n", "y4 = np.array([6.58, 5.76, 7.71, 8.84, 8.47, 7.04, 5.25, 12.50, 5.56, 7.91, 6.89])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solution" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "# Your solution" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3.3 Central limit theorem (10 points)\n", "\n", "The central limit theorem states that the sum of a large number or random variables (regardless of their distribution) tends towards a Gaussian.\n", "Consider three very different PDFs from the `scipy.stats` library: a Burr (`burr`), a double $\\Gamma$ (`dgamma`) and a $\\chi^2$ (`chi2`) functions.\n", "Below is the code to import them and initialise them with some shape parameters." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from scipy.stats import burr, chi2, dgamma\n", "pdf_chi2 = chi2.freeze(df=3)\n", "pdf_dgamma = dgamma.freeze(a=2)\n", "pdf_burr = burr(c=3, d=20, loc=-4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that you can get the value of the pdf with the method `.pdf` (e.g. `pdf_chi2.pdf(x)`) and you can generate random variates with the method `.rvs` (e.g. `pdf_chi2.rvs(N)` to generate N random variates).\n", "\n", "a) Firstly show the shape of these PDFs by plotting them in the range $[-10, 10]$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "b) Now show that the central limit theorem is valid by generating 100 random values from each of the three PDFs and summing them up. Repeat this procedure 1000 times (or more) and save the sums in an array `G`. Plot an histogram of `G`. What is the shape of its distribution? Is it the one you would expect?" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solution" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Your Solution" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3.4 Weighted mean (10 points)\n", "Consider a virtual measurement of a quantity $a$, with true value $a_{true}=3$ and the measured values being distributed around this value in a Gaussian with width $\\sigma_1=1$\n", "\n", "a) Write a program which randomly generates 1000 such measurements (draw 1000 numbers from the Gaussian). Plot the result in a histogram. Calculate the root-mean-squared deviation (rms) of the measured points.\n", "\n", "b) Consider a second measurement of the quantity with the better resolution $\\sigma_2=0.5$. For pairs of the measurements, calculate the arithmetic mean of them. Fill the arithmetic mean for 1000 pairs of such measurements in another histogram and compare the rms. Compare the rms also to a calculation. Why does taking the mean not improve the resolution?\n", "\n", "c) Now use the weighted mean presented in the lecture: $$\\hat{x}=\\frac{x_1/\\sigma^2_1+x_2/\\sigma^2_2}{1/\\sigma^2_1+1/\\sigma^2_2}$$\n", "Calculate the expected resolution of the weighted mean and compare it to your rms.\n", "\n", "Hint: You can for example use the function numpy.random.normal() to get the random numbers, check the function matplotlib.pyplot.hist() for creating histograms" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Solution" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Your Solution" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.8.3" } }, "nbformat": 4, "nbformat_minor": 4 }