mirror of
https://github.com/TeachBooks/manual.git
synced 2026-11-12 22:35:51 +00:00
Added missing files now that sphinx-download-replacer is added
This commit is contained in:
Tom van Woudenberg
parent
0c00a4cdce
commit
2474bb4899
2 changed files with 424 additions and 0 deletions
203
book/examples/example_gumbel_stripped.ipynb
Normal file
203
book/examples/example_gumbel_stripped.ipynb
Normal file
|
|
@ -0,0 +1,203 @@
|
|||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"source": [
|
||||
"## Gumbel Distribution Exercise\n",
|
||||
"\n",
|
||||
"Imagine you are concerned with the concentration of an airborne contaminant, $X$, measured in ppm. A previous benchmark study estimates that there is a 10% and 1% chance, respectively, of exceeding 4 and 10 ppm, respectively. You have been asked by the regulatory agency to estimate the contaminant concentration with 0.1% probability of being exceeded. Prior studies suggest that a Gumbel distribution can be used to model contaminant concentration, given by the CDF:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"F_X(x) = \\exp{\\Bigg( -\\exp{\\Big(-\\frac{x-\\mu}{\\beta}\\Big)}\\Bigg)}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"Using the cell blocks below as a guide (and also to check your analysis): Task 1) find the parameters of a Gumbel distribution that matches the information provided above, then, Tasks 2-3) use it to estimate the concentration with exceedance probability 0.1%.\n",
|
||||
"\n",
|
||||
"To complete this assignment, you can use `numpy`, `matplotlib` and from the `math` library, `log` and `e`"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"import matplotlib.pyplot as plt\n",
|
||||
"from math import log, e\n",
|
||||
"import numpy as np"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"**Task 0:** fill in the appropriate fitting points:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"x_1 = _\n",
|
||||
"p_1 = _\n",
|
||||
"x_2 = _\n",
|
||||
"p_2 = _"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"**Task 1:** derive the distribution parameters:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def gumbel_2_points(x_1, p_1, x_2, p_2): \n",
|
||||
" \"\"\"Compute Gumbel distribution parameters from two points of the CDF.\n",
|
||||
" Arguments:\n",
|
||||
" x_1 (float): point one\n",
|
||||
" p_1 (float): cumulative probability for point 1\n",
|
||||
" x_2 (float): point two\n",
|
||||
" p_2 (float): cumulative probability for point 2\n",
|
||||
" \"\"\"\n",
|
||||
"\n",
|
||||
" # YOUR CODE GOES HERE #\n",
|
||||
" beta = _\n",
|
||||
" mu = _\n",
|
||||
" #######################\n",
|
||||
"\n",
|
||||
" return mu, beta"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The cells below will print your parameter values and create a plot to help confirm you have the right implementation."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def plot_distribution(mu, beta,\n",
|
||||
" x_1, p_1, x_2, p_2):\n",
|
||||
" plt.title(\"Gumbel Distribution, $1-F_X(x)$\")\n",
|
||||
" plt.xlabel(\"Contaminant Concentration, $X$ [ppm]\")\n",
|
||||
" plt.ylabel(\"Exceedance Probability [--]\")\n",
|
||||
" plt.grid(color='black', linestyle='-', linewidth=0.3)\n",
|
||||
" x_axis = np.arange(0, 20, 0.1)\n",
|
||||
" plt.plot(x_axis, np.vectorize(gumbel_distribution)(x_axis), linewidth=2)\n",
|
||||
" plt.plot(x_1, 1 - p_1, 'ro')\n",
|
||||
" plt.annotate(\"Point 1\", (x_1 + 0.4, 1 - p_1 + 0.02))\n",
|
||||
" plt.plot(x_2, 1 - p_2, 'ro')\n",
|
||||
" plt.annotate(\"Point 2\", (x_2 + 0.4, 1 - p_2 + 0.001))\n",
|
||||
" plt.yscale(\"log\") \n",
|
||||
" plt.show()"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"mu, beta = gumbel_2_points(x_1, p_1, x_2, p_2)\n",
|
||||
"\n",
|
||||
"print(f\"Your mu: {mu:0.5f}\\nYour beta: {beta:0.5f}\")\n",
|
||||
"\n",
|
||||
"gumbel_distribution = lambda x: 1 - e**(-e**(-(x - mu)/beta))\n",
|
||||
"\n",
|
||||
"plot_distribution(mu, beta, x_1, p_1, x_2, p_2)\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"**Task 2:** write a function to solve for the random variable value (it will be tested for you with the Check answer button, along with the distribution parameters)."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def find_x_with_probability_p(p): \n",
|
||||
" \"\"\" Compute point in the gumbel distribution for which the CDF is p\n",
|
||||
" Use the variables mu and beta defined above!\n",
|
||||
" Hint: they have been defined globally, so you don't need to \n",
|
||||
" include them as arguments.\n",
|
||||
" \"\"\"\n",
|
||||
" # YOUR CODE GOES HERE #\n",
|
||||
" x = _\n",
|
||||
" #######################\n",
|
||||
"\n",
|
||||
" return x"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"**Task 3:** use the function `find_x_with_probability_p` to evaluate the random variable value with exceedance probabiltiy of 0.001."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# YOUR CODE GOES HERE #\n",
|
||||
"x = _\n",
|
||||
"#######################\n",
|
||||
"\n",
|
||||
"print(f\"The value of x with probability of exceedance 0.001 os {x:0.5f}\")"
|
||||
]
|
||||
}
|
||||
],
|
||||
"metadata": {
|
||||
"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.9.13"
|
||||
}
|
||||
},
|
||||
"nbformat": 4,
|
||||
"nbformat_minor": 2
|
||||
}
|
||||
221
book/examples/example_quiz_interactive_stripped.ipynb
Normal file
221
book/examples/example_quiz_interactive_stripped.ipynb
Normal file
|
|
@ -0,0 +1,221 @@
|
|||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Weighted least-squares estimation\n",
|
||||
"In ordinary least-squares estimation, we assume that all observations are equally important. In many cases this is not realistic, as observations may be obtained by different measurement systems, or under different circumstances. We want our methodology for least-squares estimation to be able to take this into account.\n",
|
||||
"\n",
|
||||
"We achieve this goal by introducing a weight matrix in the minimization problem:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\underset{\\mathrm{x}}{\\min} {\\mathrm{(y-Ax)^T W(y-Ax)}}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"In the unweighted least-squares approach, we arrive at the normal equation by pre-multiplying both sides of $\\mathrm{y=Ax}$ with the transpose of the design matrix $\\mathrm{A^T}$:\n",
|
||||
"\n",
|
||||
"$$ \n",
|
||||
"\\mathrm{y=Ax} \\; \\rightarrow \\; \\mathrm{A^T\\; y = A^T\\; A x }\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"In the weighted least-squares approach, we now need to add weight matrix $W$ to this pre-multiplication factor, i.e., $ \\mathrm{A^T W}$, to obtain the normal equation:\n",
|
||||
"\n",
|
||||
"$$ \n",
|
||||
"\\mathrm{y=Ax} \\; \\rightarrow \\; \\mathrm{A^T W \\; y = A^TW \\; A x}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"The normal matrix is now defined as $\\mathrm{N=A^T W A}$. From this, assuming that the normal matrix $N$ is invertible (non-singular) we find the weighted least-squares estimate $ \\mathrm{\\hat{x}} $, \n",
|
||||
"\n",
|
||||
"$$ \n",
|
||||
"\\begin{align*}\n",
|
||||
"\\mathrm{\\hat{x}} &= \\mathrm{(A^T W A)^{-1} A^T W y} \\\\\n",
|
||||
" &= \\arg \\underset{\\mathrm{x}}{\\min} {\\mathrm{(y-Ax)^T W(y-Ax)}}\n",
|
||||
"\\end{align*}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"We also find the derived estimate $ \\mathrm{\\hat{y}} $ and $ \\mathrm{\\hat{\\epsilon}} $:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\mathrm{\\hat{y} = A \\hat{x} = A (A^T W A )^{-1} A^T W y}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
" \\mathrm{\\hat{\\epsilon} = y - \\hat{y}= y - A \\hat{x} = y-A (A^T W A )^{-1} A^T W y = (I- A(A^T W A )^{-1} A^T W) y}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"### Quiz question\n",
|
||||
"<iframe src=\"https://tudelft.h5p.com/content/1292046737060674407/embed\" aria-label=\"WLS_1\" width=\"1088\" height=\"637\" frameborder=\"0\" allowfullscreen=\"allowfullscreen\" allow=\"autoplay *; geolocation *; microphone *; camera *; midi *; encrypted-media *\"></iframe><script src=\"https://tudelft.h5p.com/js/h5p-resizer.js\" charset=\"UTF-8\"></script>\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"### Video\n",
|
||||
"<iframe width=\"560\" height=\"315\" src=\"https://www.youtube.com/embed/iJmkkz37EuU\" title=\"YouTube video player\" frameborder=\"0\" allow=\"accelerometer; autoplay; clipboard-write; encrypted-media; gyroscope; picture-in-picture\" allowfullscreen></iframe>\n",
|
||||
"\n",
|
||||
"### Discussion on the weight matrix\n",
|
||||
"The weight matrix $\\mathrm{W}$ expresses the (relative) weights between the observations. It is always a square matrix. The size of the weight matrix depends on the number of observations, $m$. The size of the weight matrix is $m\\times m$. \n",
|
||||
"\n",
|
||||
"If it is a unit matrix ($\\mathrm{W=I}$), this implies that all observations have equal weight. Note that in this case the equations are equal to the ordinary least-squares solution. \n",
|
||||
"\n",
|
||||
"If it is a diagonal matrix, with different values on the diagonal, this implies that entries with a higher value correspond to the observations which are considered to be of more importance. If the weight matrix has non-zero elements on the off-diagonal positions, this implies that (some) observations are correlated.\n",
|
||||
"\n",
|
||||
"### Weighted least-squares estimator: properties\n",
|
||||
"\n",
|
||||
"Until now, we have looked at the weighted least-squares solution of a single *realization* of the observations, where generally we assume that it is a realization of the *random* observable vector $Y$, since measurements are affected by random errors. As such it follows the the weighted least-squares *estimator* is given by:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\hat{X} = \\mathrm{(A^T W A )^{-1} A^T W} Y\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"This estimator has two important properties: it is *linear* and *unbiased*.\n",
|
||||
"\n",
|
||||
"The linearity property is due to the fact that $\\hat{X}$ is a linear function of the observables $Y$.\n",
|
||||
"\n",
|
||||
"The unbiased property means that the expectation of $\\hat{X}$ is equal to the true (but unknown) $\\mathrm{x}$. This can be shown as follows:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\mathbb{E}(\\hat{X}) = \\mathrm{(A^T W A )^{-1} A^T W} \\mathbb{E}(Y) = \\mathrm{(A^T W A )^{-1} A^T W Ax = x}\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"This a very desirable property. It applies that if we would repeat the measurements many times to obtain a new estimate, the *average of the estimated* values would be equal to the truy values."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Exercise\n",
|
||||
"\n",
|
||||
"You have a time series of 8 measurements and fit a model assuming a linear trend (constant velocity)."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"import numpy as np\n",
|
||||
"import matplotlib.pyplot as plt"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# times of observation [months]\n",
|
||||
"t = np.arange(8)\n",
|
||||
"\n",
|
||||
"# observed heights [m]\n",
|
||||
"y = [1.39, 1.26, 1.48, 4.03, 5.89, 5.14, 6.47, 7.64]\n",
|
||||
"\n",
|
||||
"# number of observations\n",
|
||||
"m = len(t)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Fill in the correct $\\mathrm{A}$-matrix."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# design matrix\n",
|
||||
"A = np.column_stack((np.ones(m), t))"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Define the weight matrix for $\\mathrm{W}=I_m$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# Weight matrix for case 1\n",
|
||||
"W_1 = ?"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Define the weight matrix with the weight of first 3 observations is five times as large as the rest"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# Weight matrix for case 2\n",
|
||||
"W_2 = ?"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Define the weight matrix with the weight of 4th and 5th observation is five times as large as the rest"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {
|
||||
"tags": []
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# Weight matrix for case 3\n",
|
||||
"W_3 = ?"
|
||||
]
|
||||
}
|
||||
],
|
||||
"metadata": {
|
||||
"kernelspec": {
|
||||
"display_name": "base",
|
||||
"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.9.18"
|
||||
},
|
||||
"orig_nbformat": 4
|
||||
},
|
||||
"nbformat": 4,
|
||||
"nbformat_minor": 2
|
||||
}
|
||||
Loading…
Add table
Add a link
Reference in a new issue