## Exercise 1 (5 points): Independent Component Analysis [Python]


While PCA (and SVD for that matter) tries to seek new axes which are designed to explain the variance in the data and hence remove correlations between feature dimensions (the covariance matrix of PCA transformed data is equal to the identity matrix), one can go even one step further and also remove statistical dependence between the data which is stronger than any correlation
removal technique. As a result, we retrieve the independent components which produced the original data. The algorithm designed for this purpose is called independent component analysis (ICA) and in this exercise, we are going to take a look at it.

Note: this exercise is based on the introduction to ICA in [1]. In case you have any questions about the theoretical background of the calculations we are going to perform here, feel free to read the corresponding section in this (very well written) paper.

In ICA, we assume a linear mixture 𝐴 between the sources in 𝑆 producing the output 𝑋. A common example is the so-called cocktail party problem [1] where we have two independent audio sources (e.g. music and voice) and two different recordings of the sources (e.g. two microphones). A linear combination of the two produced signals reaches each microphone so that our recordings 𝑋 are mixed together $as^1$ <br/><br/>
\begin{equation}X = SA \end{equation}

where both, 𝑋 and 𝑆, are matrices with the columns representing variables (here audio sources) and the rows observations (audio samples over time). In our case, we are using two sources which were manually combined over 10k samples. Hence, $X, S \in \mathbb{R}^{10000×2}$ and $A \in \mathbb{R}^{2×2}$. You only get the two “recorded” audio signals in 𝑋 and your task is to find the original source signals in $S$.

Mathematically, we need the inverse of the matrix 𝐴 to recover 𝑆. A common approach is to
simplify this by applying SVD first on 𝐴 so that the problem reduces to find three individual
matrices [1]
    \begin{equation}S = XA^{−1} = X(U\Sigma VT)^−1 = XV\Sigma^{-1}U T.\end{equation}

In the following, plot the result after each calculation. A summary of all plots which you should have in the end is shown in Figure 1. If you want to play the audio files in your Jupyter notebook, you can use the Audio function from the IPython.display package. For the test script, please pay also attention to Table 1.

1. The two audio signals are stored in the files $mixedA.wav$ and $mixedB.wav$. Load both files into your Python program and store them in one matrix (this is 𝑋). Note: the WAV files use the 32-bit floating-point format so that the loaded signal values should already be in the range [−1; 1].


2. Calculate the centralized data $X_c$ by subtracting the (column) mean $\mu_i$

\begin{equation}X_c = X − (\mu_1 \hspace{0.2cm} \mu_2) \end{equation}

3. We approach Equation 1 from the left side. It turns out that 𝑉 contains the eigenvectors of $X_c^T X_c$ in the columns and $\Sigma$ stores the square root of the eigenvalues of $X_c^TX_c$ on the diagonal (even though SVD was applied to 𝐴 and not $X_c$). Hence, we need to apply PCA to $X_c$ first and use the result to transform the data
\begin{equation}X_w = X_c V \Sigma^{-1} \end{equation}<br/>
This process is called whitening or sphering [1] of the data and there is also an option in sklearn.decomposition.PCA which you can use to calculate $X_w$.

4. The next step is to calculate $U^T$. To do so, we exploit the fact that $U^T$ is essentially a rotation matrix. This leads to a new (and not really scalable) approach to find $U^T$: we simply try a few rotation angles 𝜃 and choose the one which fits best.
To find out if a rotation matrix leads to a good or bad result, we exploit the fact of assumed independence of the signals and use a measure which tells us how independent our data currently are. This can be achieved by minimizing the so-called multi-information [1] of the rotated data2
\begin{equation}\theta^* = argmin_{\theta} \sum_{i=1}^2H ((X_wR(\theta))_i)\end{equation}
where H calculates the entropy based on the 𝑖-th column of the data transformed via the rotation matrix $R(\theta)$.

    a) Calculating the entropy of a vector 𝒙 ∈ ℝ10000 is not really applicable since it is not unlikely that most of the values occur only once leading lead to a very high entropy value which may not even change for different 𝜃 values. One way to overcome this issue is by grouping similar values together and then calculate the entropy from these groups. More concretely, this can be achieved by calculating a normalized histogram so that we have a discrete probability 𝑃(𝑖) for each of the 𝐵 bins. Then, we can simply apply the formula of the entropy 
\begin{equation}H = \sum_{i=1}^BP(i).log_2(P(i))\end{equation}
Write a function gentropy(values) which calculates the entropy based on this histogram approach. The number of bins 𝐵 is crucial here since it determines how many values are grouped together and are hence considered equal for the entropy calculation. You can set it to 𝐵 = 500.

    b) Test 180 different angular values from the range [0; 𝜋] and calculate the multiinformation with the help of your gentropy function and extract the angular value $\theta^*$ with the lowest corresponding multi-information.

5. We now have all ingredients together to reconstruct the original signals.

    a) Apply the final transformation matrix
\begin{equation}X_r = X_wU^T = X_wR(\theta^*)\end{equation}

    b) We can’t reconstruct the original amplitude values with ICA. Hence, the reconstruction is only correct up to a scaling factor. To avoid overdrive, scale the reconstructed signals $X_r$ to the range [−1; 1] via min-max scaling.

## Exercise 2 (5 points): Autoencoder [Python]

An autoencoder is an unsupervised learning technique based on neural networks which can be used to learn new features. The general idea is to construct a network hierarchy with an equal number of neurons in the input as well as in the output layer and use the input also as the teaching signal. Further, a bottleneck is introduced in the middle of the network hierarchy so that a simple forwarding of the signal is not possible. Instead, the network is forced to learn representative features to recover the data even when only a subset of the information is available.

There are three common use cases of an autoencoder: data denoising, variation learning and dimensionality reduction. In this exercise, we are going to apply the latter to a small example in order to explore the basics of autoencoders. We have already covered some dimensionality reduction techniques and it turns out that there is a relationship between autoencoders and PCA. We take this as a starting point and then extend to the general concept of non-linear PCA (with similar results as kernel PCA), i.e. instead of lines as main axes of our new coordinate system (eigenvectors) we now may end up with more arbitrary shapes which can reflect the data distribution better.

This exercise comes with some predefined code in the Autoencoder.ipynb Jupyter notebook. Use this file as the starting point for your implementation. As you can see, we are using a two-dimensional dataset here so that our projections are going to be one-dimensional. For the test script, please pay also attention to Table 2.

1. We can construct a very simple autoencoder for our dataset by using two input neurons, one neuron in the coding layer (bottleneck) and two output neurons as well as linear transfer functions. Then, we basically compute PCA and the output of the coding layer corresponds to the projections (also called codings) and the output of the network to the recovered data points based on these projections. 

    a) This type of autoencoder is already implemented. Familiarize yourself with the code and execute it.
    
    b) As a comparison, compute PCA also in the normal way and extract the first eigenvector $\mu_1$.
    
    c) Visualize the data, the output of the autoencoder, the line of the eigenvector $\mu_1$ and some contour lines which show how the data are projected onto the new axis. That is, all points on one contour line map to the same projection value. A possible result is shown in Figure 3 (a)

2. We now want to change the autoencoder so that it performs non-linear PCA instead.

    a) Implement the hierarchy as shown in Figure 4 with the following notes in mind:
        • The en- and decoder now include an additional hidden layer with 20 neurons each.
        • The hidden layers use tanh(𝑥) as transfer functions and the code as well as the output layer use linear transfer
        functions.
        • It is recommended to constrain the autoencoder network in some way or another to get good codings. Here, we want to regularize all weights by adding an L1 penalty term to the error function
\begin{equation}E = E_0 + \lambda \sum_w |w| \end{equation}
        with the base loss $E_0$ (MSE) and the regularization term 𝜆 = 0.002. Note: you need to find out yourself how this can be implemented in TensorFlow.
        
    b) Visualize the result in the same way as before (without the eigenvector $\mu_1$). A possible solution is shown in Figure 3 (b)

3. Do you think that it would be a good idea to use a non-linear transfer function (like tanh(𝑥)) in the output layer (for our dataset)? If not, why?

4. Figure 2 shows a histogram of the projected values (the codings) of both autoencoders. As you can see, the histogram of the non-linear autoencoder is more balanced whereas in the linear case we have relatively clear peaks. Explain why we observe these differences.

775b8724e08eba7f0321b9d3989f138fb06d713e
ipynb
Jupyter Notebook
11_ica_autoencoder/Autoencoder.ipynb
jhinga-la-la/pattern-recognition-course
7ad4f70b2c427f3c37f59f47768b90371873823c
Apache-2.0
11_ica_autoencoder/Autoencoder.ipynb
jhinga-la-la/pattern-recognition-course
7ad4f70b2c427f3c37f59f47768b90371873823c
Apache-2.0
11_ica_autoencoder/Autoencoder.ipynb
jhinga-la-la/pattern-recognition-course
7ad4f70b2c427f3c37f59f47768b90371873823c
Apache-2.0
Qwen/Qwen-72B
1. YES
2. YES