|
1 | 1 | # Sample Statistics and Gaussians Exercise |
2 | 2 |
|
3 | | -### Exercise 1: Mean and Variance |
4 | | -This exercise first explores sample statistics like mean and variance. |
| 3 | +Today's exercise first explores sample statistics like mean and variance. |
5 | 4 | It continues with the definition of the gaussian distributions. |
6 | 5 | It ends with applying gaussian mixture models (gmm) and the expectation-maximization algorithm used to optimize gmm. |
7 | 6 |
|
| 7 | +### ⊙ Task 1: Mean and Variance |
| 8 | + |
8 | 9 | - To get started, take a look at the `src/sample_mean_corr_rhein.py` file. |
9 | 10 |
|
10 | 11 | Implement functions to compute sample mean and the standard deviation. |
11 | | -Use |
12 | | - |
13 | | -$$ \hat{\mu} = \frac{1}{n} \sum_{i=1}^n x_i , $$ |
| 12 | +1. Use |
14 | 13 |
|
15 | | -to calculate the mean and |
| 14 | + $$ \hat{\mu} = \frac{1}{n} \sum_{i=1}^n x_i , $$ |
16 | 15 |
|
17 | | -$$ \hat{\sigma} = \sqrt{\frac{1}{n-1} \sum_{i=1}^n (x_i - \hat{\mu})^2} $$ |
| 16 | + to calculate the mean. |
18 | 17 |
|
19 | | -to compute the standard deviation. $x \in \mathbb{R}$ denotes individual sample elements, and $n \in \mathbb{N}$ the size of the sample. |
| 18 | +2. Use |
20 | 19 |
|
21 | | -Return to the Rhine data-set. Load the data from `./data/pegel.tab`. Compute the water level mean and standard deviation before and after the year 2000. |
| 20 | + $$ \hat{\sigma} = \sqrt{\frac{1}{n-1} \sum_{i=1}^n (x_i - \hat{\mu})^2} $$ |
22 | 21 |
|
| 22 | + to compute the standard deviation. |
| 23 | +$x_i \in \mathbb{R}$ for $i \in \{1, ... , n\}$ denotes individual sample elements, and $n \in \mathbb{N}$ the size of the sample. |
| 24 | +Don't use the pre-build functions np.mean() or np.std() to solve these tasks. |
23 | 25 |
|
| 26 | +Return to the Rhine data-set and go to the `__main__`-function. |
| 27 | +The data from `./data/pegel.tab` is already loaded and processed and is ready to be used. |
| 28 | +3. Compute the water level mean and standard deviation before the year 2000. |
| 29 | +4. And now compute the water level mean and standard deviation after the year 2000. |
24 | 30 |
|
25 | | -### Exercise 2: Autocorrelation |
| 31 | +### ⊙ Task 2: Autocorrelation |
26 | 32 | We now want to use autocorrelation to analyse the discrete time signal of the rhine level measurements. Implement the `auto_corr` function in `src/sample_mean_corr_rhein.py`. It should implement the engineering version without the normalization and return the autocorrelation |
27 | 33 |
|
28 | 34 | $$ R_{xx} = (c_{-N+1},\ldots,c_{1}, c_0, c_{1}, \ldots, c_{N-1}) $$ |
|
31 | 37 |
|
32 | 38 | $$ c_{k} = \sum_{t=1}^{N-|k|} n_t n_{t + |k|}$$ |
33 | 39 |
|
34 | | -with $n$ the normalized version of your signal of length $N$. The time shift $k$ moves from $-(N-1)$ to $N-1$. Therefore, the resulting array has a length of $2N-1$. For example the autocorrelation of an input signal $x=(2,3,-1)$ is $R_{xx}=(c_{-2}, c_{-1}, c_0, c_1, c_2)=(-2, 3, 14, 3, -2)$ and is symmetrical. Make sure that you normalize the signal *before* giving it to your `auto_corr` function. Once you have checked your implementation using `nox -s test`, you can use `np.correlate` for efficiency. Plot the autocorrelation for the rhine level measurements since 2000. |
35 | | -Normalize your data via |
| 40 | +with $n$ the normalized version of your signal of length $N$. The time shift $k$ moves from $-(N-1)$ to $N-1$. Therefore, the resulting array has a length of $2N-1$. |
| 41 | + |
| 42 | +For example the autocorrelation of an input signal $x=(2,3,-1)$ is $R_{xx}=(c_{-2}, c_{-1}, c_0, c_1, c_2)=(-2, 3, 14, 3, -2)$ and is symmetrical. |
| 43 | + |
| 44 | +>>> In the following table you can see an illustrative depicition on how the $c_k$'s are calculated. |
| 45 | +The header contains the input signal x padded with 0's on its sides. |
| 46 | +For autocorrelation, we compute correlation between $x$ and $x$ itself. |
| 47 | +So in visual terms, we slide $x$ from left to right across itself. |
| 48 | +At each step we compute one $c_k$ by first multiplying the numbers that are aligned with the input signal in the header. Then, these products will be summed up. The result is written in the respective cell of the last column. |
| 49 | + |
| 50 | +| 0 | 0 | 2 | 3 | -1 | 0 | 0 | $c_k$ | |
| 51 | +| -------- | ------- | ------- | ------- | ------- | ------- | ------- | ------- | |
| 52 | +| 2 | 3 | -1 | | | | | 0 + 0 - 2 = -2 |
| 53 | +| | 2 | 3 | -1 | | | | 0 + 6 - 3 = 3 |
| 54 | +| | | 2 | 3 | -1 | | | 4 + 9 + 1 = 14 |
| 55 | +| | | | 2 | 3 | -1 | | 6 - 3 + 0 = 3 |
| 56 | +| | | | | 2 | 3 | -1 | -2 + 0 + 0 = -2 |
| 57 | +>>>As you can see, when reading from top to bottom we get the correct solution $R_{xx}=(-2, 3, 14, 3, -2)$. |
| 58 | +
|
| 59 | +So here are your tasks: |
| 60 | +1. Implement the `auto_corr` function as described above. |
| 61 | + >>> The function expects $x$ to be normalized. That means, that no normalization is done inside the `auto_corr`-function. Instead you normalize the input signal before calling `auto_corr`. |
| 62 | +2. Check your implementation using `nox -s test`. |
| 63 | + If the test passes you can use `np.correlate` for efficiency in the following exercises! |
36 | 64 |
|
37 | | -$$ n_{t} = \frac{x_{t} - \hat{\mu}}{\hat{\sigma}} ,$$ |
| 65 | +Now go back to the `__main__`-function and consider the Rhine data-set after 2000. |
| 66 | +3. Normalize the data of the Rhine level measurements since 2000 via |
38 | 67 |
|
39 | | -for all $t$ measurements until the signal length N. Before running the autocorrelation computation. Compare the autocorrelation to a random signal from `np.random.randn` by plotting both results with `plt.plot`. |
| 68 | + $$ n_{t} = \frac{x_{t} - \hat{\mu}}{\hat{\sigma}} ,$$ |
40 | 69 |
|
| 70 | + for all $t$ measurements until the signal length N. |
| 71 | +4. Compute and plot the autocorrelation for the Rhine level measurements since 2000. |
41 | 72 |
|
42 | | -### Exercise 3: Distributions |
| 73 | +Now we want to compare this autocorrelation to the one of a random signal. |
| 74 | +5. Create a random signal from `np.random.randn` that has the same shape as the Rhine level measurements |
| 75 | +6. Normalize the random signal. |
| 76 | +7. Compute the autocorrelation of the normalized, random signal. |
| 77 | +8. Plot both autocorrelations of the Rhine level measurements and the random signal using `plt.plot` and compare the results. |
43 | 78 |
|
44 | 79 |
|
45 | | -- Consider the `src/plot_gaussian.py` module. Implement the `gaussian_pdf` function. |
| 80 | +### ⊙ Task 3: Distributions |
46 | 81 |
|
47 | | -In one dimension gaussian probability density function is defined as |
| 82 | +1. Consider the `src/plot_gaussian.py` module. Implement the `gaussian_pdf` function. |
48 | 83 |
|
49 | | -$$\phi_1(x | \mu, \sigma) = \frac{1}{\sqrt{2\pi \sigma^2}} \exp({-\frac{1}{2}(\frac{x - \mu}{\sigma})^2}) .$$ |
| 84 | + In one dimension gaussian probability density function is defined as |
50 | 85 |
|
51 | | -$\pi \in \mathbb{R}$ denotes Pi, $\mu \in \mathbb{R}$ the mean and $\sigma \in \mathbb{R}$ the standard deviation for a random variable $X$. $e^x$ denotes the exponential function. $X$ having a gaussian pdf is described as gaussion or normal distribution $\mathcal{N}$. Explore the behavior of $\mathcal{N}(\mu, \sigma)$ for different values of $\mu$ and $\sigma$. |
| 86 | + $$\phi_1(x | \mu, \sigma) = \frac{1}{\sqrt{2\pi \sigma^2}} \exp({-\frac{1}{2}(\frac{x - \mu}{\sigma})^2}) .$$ |
52 | 87 |
|
| 88 | + $\pi \in \mathbb{R}$ denotes Pi, $\mu \in \mathbb{R}$ the mean and $\sigma \in \mathbb{R}$ the standard deviation for a random variable $X$. |
| 89 | + $e^x$ denotes the exponential function. |
| 90 | + A random variable $X$ having a gaussian pdf is described as gaussion or normal distribution $\mathcal{N}$. |
| 91 | + >>> Remark: In the notation $\phi_1(x | \mu, \sigma)$, x is the variable that is plugged into the function and $\mu$ and $\sigma$ are parameters which are needed to define the function and that are determined beforehand. |
| 92 | +2. Explore the behavior of $\mathcal{N}(\mu, \sigma)$ for different values of $\mu$ and $\sigma$. |
| 93 | + The Code for plotting the pdf's is already given. |
53 | 94 |
|
54 | | -- Consider the `src/mixture_concpets.py` module. |
55 | | -Implement a two-dimensional gaussian pdf following, |
56 | 95 |
|
57 | | -$$ \phi_2(\mathbf{x} | \mu_g, \Sigma_g) = \frac{1}{\sqrt{(2\pi)^2 \| \Sigma_g \|}} \exp({-\frac{1}{2}(\mathbf{x}-\mu_g)^T \Sigma_g^{-1}(\mathbf{x}-\mu_g)}).$$ |
58 | 96 |
|
59 | | -$\mu_g \in \mathbb{R}^2$ denotes the two dimensional mean vector, $\Sigma_g \in \mathbb{R}^{2\times2}$ the covariance matrix, $^{-1}$ the matrix inverse, $T$ the transpose and $g \in \mathbb{N}$ the number of the distrubtion, which will be important later. |
60 | | -Plot a 2d-bell curve with $\mu_1 = [-1.5, 2]$ and $\Sigma_1 = [[1, 0], [0, 1]]$ using the `plt.imshow` function. `np.linspace` and `np.meshgrid` will help you. |
| 97 | +3. Consider the `src/mixture_concepts.py` module. |
| 98 | +Go to the `twod_gaussian_pdf`-function and implement a two-dimensional gaussian pdf following, |
61 | 99 |
|
| 100 | + $$ \phi_2(\mathbf{x} | \mu_g, \Sigma_g) = \frac{1}{\sqrt{(2\pi)^2 \| \Sigma_g \|}} \exp({-\frac{1}{2}(\mathbf{x}-\mu_g)^T \Sigma_g^{-1}(\mathbf{x}-\mu_g)}).$$ |
62 | 101 |
|
63 | | -### Exercise 4: Gaussian mixture models (optional) |
| 102 | + $\mu_g \in \mathbb{R}^2$ denotes the two dimensional mean vector, $\Sigma_g \in \mathbb{R}^{2\times2}$ the covariance matrix, $^{-1}$ the matrix inverse, $T$ the transpose, $\| \|$ the determinant and $g \in \mathbb{N}$ the number of the distrubtion, which will be important later. |
| 103 | + |
| 104 | + - As you can see, the x-parameter of the function is a grid of shape (grid_height, grid_width, 2). |
| 105 | + That means, that you get not only one but grid_height*grid_width many 2-dimensional values that should be evaluated. |
| 106 | + It's up to you how you want to approach this task. |
| 107 | + But Broadcasting might be an elegant way to evaluate all these values at the same time (https://numpy.org/doc/stable/user/basics.broadcasting.html). In this case, you might take a look at `np.swapaxes` to deal with the transponation. |
| 108 | + - At the very end of this document, we included a diagram that depicts how the shapes should develop when using Broadcasting. This is purely optional and just for when you need some guidance regarding the relevant shapes. |
| 109 | + |
| 110 | +4. Plot a 2d-bell curve with $\mu_1 = [-1.5, 2]$ and $\Sigma_1 = [[1, 0], [0, 1]]$ using the `plt.imshow` function. `np.linspace` and `np.meshgrid` will help you. |
| 111 | + |
| 112 | + |
| 113 | + |
| 114 | +### ✪ Task 4: Gaussian mixture models (optional) |
64 | 115 |
|
65 | 116 | We can use bell-curve sums for classification! A Gaussian mixture model has the density |
66 | 117 |
|
@@ -113,3 +164,8 @@ Train a gmm to find the diabetic patients. |
113 | 164 |
|
114 | 165 |
|
115 | 166 | - Standard packages like sci-kit-learn implement GMMs, too. Take a minute to read https://scikit-learn.org/stable/modules/mixture.html . |
| 167 | + |
| 168 | +--- |
| 169 | +If you need some inspiration for broadcasting in Task 3.3: |
| 170 | + |
| 171 | + |
0 commit comments