<jupyter_start><jupyter_text>Lab 13- PCA- TSNE  PCA (Principal Component Analysis)Analiza głównych składników jest szybką metodą redukcji wymiarowości danych.Polega na znalezieniu liniowej transformacji zbioru zmiennych początkowych w mniej liczny zbiór zmiennych, zwanych składowymi głównymi. Jest to pewien sposób kompresji danych, dlatego metodę tę można wykorzystać w wielu dziedzinach nauki.**Algorytm PCA:**1. Standaryzacja danych2. Wyznaczenie macierzy kowariancji $Σ$ między zmiennymi początkowymi3. Wyznaczenie wartości własnych macierzy kowariancji4. Wybranie $k$ największych wartości własnych i wyznaczenie dla nich wektorów własnych5. Utworzenie macierzy przekształcenia liniowego $W$, bazującej na wektorach własnych6. Przekształcenie zmiennych początkowych według wzoru $Y = WX$![](./fig/pca_matrix.png)![](./fig/pca_gif.gif)  Przykład dla sztucznych danych<jupyter_code>%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns; sns.set()
rng = np.random.RandomState(1)
X = np.dot(rng.rand(2, 2), rng.randn(2, 200)).T
X.shape
plt.figure(figsize=(8,5))
plt.scatter(X[:, 0], X[:, 1])
plt.show()<jupyter_output><empty_output><jupyter_text>Jeden komponent<jupyter_code>from sklearn.decomposition import PCA
pca = PCA(n_components=1)
pca.fit(X)
# Kierunki zmiennych PCA
pca.components_
# Wyjaśniona wariancja
print(f'Total variance: {np.sqrt(np.std(X)):.4f}')
pca.explained_variance_
# Procent wyjaśnionej wariancji
pca.explained_variance_ratio_
def draw_vector(v0, v1, ax=None):
    ax = ax or plt.gca()
    arrowprops=dict(arrowstyle='->',
                    linewidth=4,
                    shrinkA=0, shrinkB=0, color='black')
    ax.annotate('', v1, v0, arrowprops=arrowprops)

# plot data
plt.figure(figsize=(8,5))
plt.scatter(X[:, 0], X[:, 1])
for length, vector in zip(pca.explained_variance_, pca.components_):
    v = vector * 3 * np.sqrt(length)
    draw_vector(pca.mean_, pca.mean_ + v)
plt.show()<jupyter_output><empty_output><jupyter_text>2 komponenty<jupyter_code>from sklearn.decomposition import PCA
pca = PCA(n_components=2)
pca.fit(X)
# Kierunki zmiennych PCA
pca.components_
# Wyjaśniona wariancja
print(f'Total variance: {np.sqrt(np.std(X)):.4f}')
pca.explained_variance_
# Procent wyjaśnionej wariancji
pca.explained_variance_ratio_
def draw_vector(v0, v1, ax=None):
    ax = ax or plt.gca()
    arrowprops=dict(arrowstyle='->',
                    linewidth=4,
                    shrinkA=0, shrinkB=0, color='black')
    ax.annotate('', v1, v0, arrowprops=arrowprops)

# plot data
plt.figure(figsize=(8,5))
plt.scatter(X[:, 0], X[:, 1])
for length, vector in zip(pca.explained_variance_, pca.components_):
    v = vector * 3 * np.sqrt(length)
    draw_vector(pca.mean_, pca.mean_ + v)
plt.show()<jupyter_output><empty_output><jupyter_text>Dobieranie odpowiedniej liczby komponentów<jupyter_code>import pandas as pd
hitters = pd.read_csv('hitters.csv')
hitters.head()<jupyter_output><empty_output><jupyter_text>Wydzielenie zmiennej celu i podział na zbiór treningowy i testowy<jupyter_code>from sklearn.model_selection import train_test_split

y = hitters['NewLeague']
X = hitters.drop(['Name','NewLeague','League','Division'], axis=1)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.10, random_state=0)
from sklearn.decomposition import PCA
pca = PCA().fit(X_train)

plt.figure(figsize=(9,6))
plt.plot(np.cumsum(pca.explained_variance_ratio_))
plt.xlabel('number of components')
plt.ylabel('cumulative explained variance');<jupyter_output><empty_output><jupyter_text>Zobaczmy jak poradzi sobie klasyfikacja dla surowych danych i dla PCA<jupyter_code>from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import accuracy_score
model = RandomForestClassifier(random_state=0)

# Surowe dane
y_hat = model.fit(X_train,y_train).predict(X_test)
print(f'Accuracy bez PCA: {accuracy_score(y_test, y_hat):.4f}')

# PCA - odczytane z wykresu
acc = []
for i in range(1,17):
    pca = PCA(n_components=i)
    pca.fit(X_train)
    y_hat = model.fit(pca.transform(X_train),y_train).predict(pca.fit_transform(pca.transform(X_test)))
    acc.append(accuracy_score(y_test, y_hat))
    
plt.figure(figsize=(9,6))
plt.plot(acc)
plt.xlabel('number of components')
plt.ylabel('Accuracy');<jupyter_output>Accuracy bez PCA: 0.5185<jupyter_text>TSNE Pytanie: Jak zwizualizować dane wielowymiarowe? **t-distributed Stochastic Neighbor Embedding** Cel: Pokazać ukryte zależności wielowymiarowych danych w 2D lub 3D.Idea: Chcemy, aby obserwacje podobne do siebie w wielu wymiarach były blisko siebie w podprzestrzeni.Wyróżnijmy:  **obserwacja x** - wielowymiarowy wektor cech o wymiarze zgodnym z wymiarem danych.  **mapowanie y** - dwu- lub trzy- wymiarowy wektor określający pozycje obserwacji na mapie.1. Obliczamy odległosci pomiędzy obserwacjami poprzez wyliczenie odpowiednich prawdopodobieństw warunkowych.Im bardziej obserwacje są podobne do siebie, tym większe prawdopodobieństwo. $$p_{j|i} = \frac{\exp{(-d(\boldsymbol{x}_i, \boldsymbol{x}_j) / (2 \sigma_i^2)})}{\sum_{i \neq k} \exp{(-d(\boldsymbol{x}_i, \boldsymbol{x}_k) / (2 \sigma_i^2)})}, \quad p_{i|i} = 0,$$Powyższy wzorek określa podobieństwo dwóch obserwacji z wykorzystaniem rozkładu normalnego scentrowanego w $x_{i}$. $\sigma_{i}$ wyznaczana jest osobno dla każdej cechy w zależności od jej wariancji. W `sklearn` możemy wpływać na "spłaszczenie" lub "zwężanie" rozkładu poprzez parametr `perplexity`.Aby skorzystać z symetrii wprowadzamy:  $$p_{ij} = \frac{p_{j|i} + p_{i|j}}{2N}.$$W ten sposób otrzymujemy macierz podobieństwa **P** (ang. *similarity matrix*)2. Następnie rozważamy macierz **Q** stanowiącą podobne odwzorowanie co powyżej, ale dla mapowań. Z przyczyn, które na razie pomińmy prawdopodobieństwa nie są liczone z wykorzystaniem rozkładu normalnego, ale t-Studenta o jednym stopniu swobody (rozkład Cauchy'ego):$$q_{ij} = \frac{(1 + ||\boldsymbol{y}_i - \boldsymbol{y}_j||^2)^{-1}}{\sum_{k \neq l} (1 + ||\boldsymbol{y}_k - \boldsymbol{y}_l||^2)^{-1}},$$3. Zaszliśmy daleko. Rozważamy teraz takie "umiejscowienie" obserwacji w przestrzeni o 2 lub 3 wymiarach, aby jak najbardziej zminimalizować różnicę pomiędzy dwiema macierzami (P i Q). Zmieniamy oczywiście macierz mapowań **Q**. Odpowiada to minimalizacji dywergencji Kullbacka-Leibera (https://pl.wikipedia.org/wiki/Dywergencja_Kullbacka-Leiblera, https://en.wikipedia.org/wiki/Kullback%E2%80%93Leibler_divergence):$$KL(P|Q) = \sum_{i \neq j} p_{ij} \log \frac{p_{ij}}{q_{ij}}$$  Przykład 1   Przykład 2 1. Użyjemy zbioru :https://scikit-learn.org/stable/modules/generated/sklearn.datasets.load_digits.htmlsklearn.datasets.load_digits<jupyter_code>from sklearn.manifold import TSNE
from sklearn.datasets import load_digits
from sklearn.preprocessing import scale

import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns<jupyter_output><empty_output><jupyter_text>2. Zawiera on 1797 obrazków cyfr w formacie 8x8 pikseli.<jupyter_code>digits = load_digits()
digits.data.shape<jupyter_output><empty_output><jupyter_text>3. Przykładowe obserwacje:<jupyter_code>nrows, ncols = 2, 5
plt.figure(figsize=(10,5))
plt.gray()
for i in range(ncols * nrows):
    ax = plt.subplot(nrows, ncols, i + 1)
    ax.matshow(digits.images[i,...])
    plt.xticks([]); plt.yticks([])
    plt.title(digits.target[i])
plt.show()<jupyter_output><empty_output><jupyter_text>4. Uporządkowanie danych (tylko, aby pomóc nam w wizualizacji, algorytm tego nie wymaga)<jupyter_code>X = np.vstack([digits.data[digits.target==i]
               for i in range(10)])
y = np.hstack([digits.target[digits.target==i]
               for i in range(10)])<jupyter_output><empty_output><jupyter_text>5. Wywołanie PCA i wizualizacja<jupyter_code>def scatter(x, colors):
    palette = np.array(sns.color_palette("hls", 10))

    f = plt.figure(figsize=(8, 8))
    ax = plt.subplot(aspect='equal')
    sc = ax.scatter(x[:,0], x[:,1], lw=0, s=40,
                    c=palette[colors.astype(np.int)])
    plt.xlim(-25, 25)
    plt.ylim(-25, 25)
    ax.axis('off')
    ax.axis('tight')

    # add labels
    txts = []
    for i in range(10):
        # Position of each label.
        xtext, ytext = np.median(x[colors == i, :], axis=0)
        txt = ax.text(xtext, ytext, str(i), fontsize=24)
        txts.append(txt)

    return f, ax, sc, txts

pca = PCA(n_components = 2)
digits_proj_pca = pca.fit_transform(X)
scatter(digits_proj_pca, y)
plt.show()<jupyter_output><empty_output><jupyter_text>6. Wywołanie t-SNE<jupyter_code>%%time
random_state = 1500100900
tSNE = TSNE(random_state=random_state, verbose=1)
digits_proj = tSNE.fit_transform(X)
scatter(digits_proj, y)
plt.show()
tSNE.fit(X).get_params()<jupyter_output>[t-SNE] Computing 91 nearest neighbors...
[t-SNE] Indexed 1797 samples in 0.008s...
[t-SNE] Computed neighbors for 1797 samples in 0.379s...
[t-SNE] Computed conditional probabilities for sample 1000 / 1797
[t-SNE] Computed conditional probabilities for sample 1797 / 1797
[t-SNE] Mean sigma: 8.121143
[t-SNE] KL divergence after 250 iterations with early exaggeration: 61.866585
[t-SNE] KL divergence after 1000 iterations: 0.750101<jupyter_text>7. Jak widzieliśmy powyżej, t-SNE jest kosztowne obliczeniowo nawet dla tak niewielkiego zbioru danych. Aby skrócić czas obliczeń, jeśli liczba cech jest znacząca, można połączyć oba podejścia: PCA + tSNE.<jupyter_code>pca =  PCA().fit(X)

plt.figure(figsize=(9,6))
plt.plot(np.cumsum(pca.explained_variance_ratio_))
plt.xlabel('number of components')
plt.ylabel('cumulative explained variance');<jupyter_output><empty_output><jupyter_text>8. Spróbujmy z `n_components=40`<jupyter_code>np.cumsum(PCA(n_components=40).fit(X).explained_variance_ratio_)[-1]
%%time
X_pca =  PCA(n_components=40).fit_transform(X)
random_state = 1500100900
digits_proj = TSNE(random_state=random_state).fit_transform(X_pca)
scatter(digits_proj, y)
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12658610271903323, 'fraction_numerical': 0.023464249748237662, 'mean_word_length': 5.594289508632138, 'pattern_counts': {'":': 0, '<': 55, '<?xml version=': 0, '>': 57, 'https://': 3, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '14', 'dir': 'jupyter-structured-clean-dedup', 'id': '51078', 'license': 'Apache-2.0', 'n_tokens_mistral': 4249, 'n_tokens_neox': 4126, 'n_words': 830, 'path': 'Materialy/Grupa3/lab_13/lab13.ipynb', 'repo_name': 'rrolkiewicz/2020L-WUM', 'size': '9930'}