<jupyter_start><jupyter_text>Chapter 5 - Compressing Data via Dimensionality Reduction  Overview - [Principal Component Analysis (PCA)](Unsupervised-dimensionality-reduction-via-principal-component-analysis)  - [Total and explained variance](Total-and-explained-variance)  - [Feature transformation](Feature-transformation)  - [Principal component analysis in scikit-learn](Principal-component-analysis-in-scikit-learn)<jupyter_code>from IPython.display import Image
%matplotlib inline<jupyter_output><empty_output><jupyter_text>Principal Component Analysis (PCA)- Unsupervised linear transformation, mostly used for:    - dimensionality reduction    - data de-noising    - data visualization<jupyter_code>Image(filename='./images/05_01.png', width=400)<jupyter_output><empty_output><jupyter_text>- PCA aims to find the **orthogonal** directions of **maximum variance** in the feature speace- PCA projects the original data points onto a new subspace with equal of fewer dimensions than the original space   - $x_1$ and $x_2$: original feature axes   - PC1 and PC2: principal components $$      z^{(i)} = v_1 x^{(i)}_1 + \dots + v_n x^{(i)}_n$$$$      z = Xv \in \mathbb{R}^m$$Variance:$$      \max_{v \in \mathbb{R}^n} \;\; \frac{1}{m} \sum_{i=1}^m (z^{(i)} - \mu_z)^2$$$$      \max_{v \in \mathbb{R}^n} \;\; \frac{1}{m} z^Tz$$$$      \max_{v \in \mathbb{R}^n} \;\; \frac{1}{m} (Xv)^T(Xv) = \frac{1}{m} v^T X^T X v = v^T \Sigma v$$  Rayleigh QuotientFor a symmetric matrix $\Sigma$, the Rayleigh quotient is defined as$$  R(\Sigma,v) := \frac{v^T \Sigma v}{v^T v}$$For the smallest and largest eigenvalues $\lambda_\min$ and $\lambda_\max$ of $\Sigma$, we can show that$$   \lambda_\min \le R(\Sigma, v) \le \lambda_\max$$where the minimum and the maximum are achieved by the eigenvectors corresponding to $\lambda_\min$ and $\lambda_\max$.  PC1Using the fact about the Rayleigh quotient, we will set$$    z = Xv_\max, $$where $v_\max$ is the eigenvector associated with $\lambda_\max$. Here, $v_\max$ is the first principal component (PC1), and $z$ contains the locations of $m$ points on the PC1.  PC2, PC3, ... PCkPC2 should be the direction of the maximum variance, excluding the variance already explained by the PC1. Further, PC2 should be orthogonal to the PC1.It turns out, that PC2, PC3, ..., PCk corresponds to the eigenvectors associated with the second largest, the third largest, ..., the $k$-th largest eigenvalues of $\Sigma$.<jupyter_code>import pandas as pd

df_wine = pd.read_csv('https://archive.ics.uci.edu/ml/'
                      'machine-learning-databases/wine/wine.data',
                      header=None)

df_wine.columns = ['Class label', 'Alcohol', 'Malic acid', 'Ash',
                   'Alcalinity of ash', 'Magnesium', 'Total phenols',
                   'Flavanoids', 'Nonflavanoid phenols', 'Proanthocyanins',
                   'Color intensity', 'Hue',
                   'OD280/OD315 of diluted wines', 'Proline']

df_wine.head()<jupyter_output><empty_output><jupyter_text>Splitting the data into 70% training and 30% test subsets.<jupyter_code>from sklearn.model_selection import train_test_split

X, y = df_wine.iloc[:, 1:].values, df_wine.iloc[:, 0].values

X_train, X_test, y_train, y_test = \
    train_test_split(X, y, test_size=0.3, random_state=0)<jupyter_output><empty_output><jupyter_text>Standardizing the data.<jupyter_code>from sklearn.preprocessing import StandardScaler

sc = StandardScaler()
X_train_std = sc.fit_transform(X_train)
X_test_std = sc.transform(X_test)<jupyter_output><empty_output><jupyter_text>Eigendecomposition of the covariance matrix.<jupyter_code>import numpy as np
cov_mat = np.cov(X_train_std.T)
#eigen_vals, eigen_vecs = np.linalg.eig(cov_mat)
eigen_vals, eigen_vecs = np.linalg.eigh(cov_mat)

print('\nEigenvalues \n%s' % eigen_vals)
print('Dim of X_train_std', np.shape(X_train_std))
print('Dim of cov_mat = ', np.shape(cov_mat))<jupyter_output>Dim of X_train_std (124, 13)
Dim of cov_mat =  (13, 13)<jupyter_text>**Note**: Above, I used the [`numpy.linalg.eig`](http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eig.html) function to decompose the symmetric covariance matrix into its eigenvalues and eigenvectors.    >>> eigen_vals, eigen_vecs = np.linalg.eig(cov_mat)    This is probably suboptimal. It would be better to use [`numpy.linalg.eigh`](http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.eigh.html) in such cases, which has been designed for [Hermetian matrices](https://en.wikipedia.org/wiki/Hermitian_matrix). The latter always returns real  eigenvalues; whereas the numerically less stable `np.linalg.eig` can decompose nonsymmetric square matrices, you may find that it returns complex eigenvalues in certain cases.   Total and explained variance<jupyter_code>tot = sum(eigen_vals)
var_exp = [(i / tot) for i in sorted(eigen_vals, reverse=True)]
cum_var_exp = np.cumsum(var_exp)
import matplotlib.pyplot as plt


plt.bar(range(1, 14), var_exp, alpha=0.5, align='center',
        label='individual explained variance')
plt.step(range(1, 14), cum_var_exp, where='mid',
         label='cumulative explained variance')
plt.ylabel('Explained variance ratio')
plt.xlabel('Principal components')
plt.legend(loc='best')
plt.tight_layout()
# plt.savefig('./figures/pca1.png', dpi=300)
plt.show()<jupyter_output><empty_output><jupyter_text>Feature transformation<jupyter_code># Make a list of (eigenvalue, eigenvector) tuples
eigen_pairs = [(np.abs(eigen_vals[i]), eigen_vecs[:, i])
               for i in range(len(eigen_vals))]

# Sort the (eigenvalue, eigenvector) tuples from high to low
eigen_pairs.sort(key=lambda k: k[0], reverse=True)
w = np.hstack((eigen_pairs[0][1][:, np.newaxis],
               eigen_pairs[1][1][:, np.newaxis]))
print('Matrix W:\n', w)
X_train_pca = X_train_std.dot(w)
colors = ['r', 'b', 'g']
markers = ['s', 'x', 'o']

for l, c, m in zip(np.unique(y_train), colors, markers):
    plt.scatter(X_train_pca[y_train == l, 0], 
                X_train_pca[y_train == l, 1], 
                c=c, label=l, marker=m)

plt.xlabel('PC 1')
plt.ylabel('PC 2')
plt.legend(loc='lower left')
plt.tight_layout()
# plt.savefig('./figures/pca2.png', dpi=300)
plt.show()
X_train_std[0]
X_train_std[0].dot(w)<jupyter_output><empty_output><jupyter_text>Principal component analysis in scikit-learn<jupyter_code>from sklearn.decomposition import PCA

pca = PCA()
X_train_pca = pca.fit_transform(X_train_std)
pca.explained_variance_ratio_
plt.bar(range(1, 14), pca.explained_variance_ratio_, alpha=0.5, align='center')
plt.step(range(1, 14), np.cumsum(pca.explained_variance_ratio_), where='mid')
plt.ylabel('Explained variance ratio')
plt.xlabel('Principal components')
plt.show()
pca = PCA(n_components=2)
X_train_pca = pca.fit_transform(X_train_std)
X_test_pca = pca.transform(X_test_std)
plt.scatter(X_train_pca[:, 0], X_train_pca[:, 1])
plt.xlabel('PC 1')
plt.ylabel('PC 2')
plt.show()
from matplotlib.colors import ListedColormap

def plot_decision_regions(X, y, classifier, resolution=0.02):

    # setup marker generator and color map
    markers = ('s', 'x', 'o', '^', 'v')
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')
    cmap = ListedColormap(colors[:len(np.unique(y))])

    # plot the decision surface
    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution),
                           np.arange(x2_min, x2_max, resolution))
    Z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)
    Z = Z.reshape(xx1.shape)
    plt.contourf(xx1, xx2, Z, alpha=0.4, cmap=cmap)
    plt.xlim(xx1.min(), xx1.max())
    plt.ylim(xx2.min(), xx2.max())

    # plot class samples
    for idx, cl in enumerate(np.unique(y)):
        plt.scatter(x=X[y == cl, 0], 
                    y=X[y == cl, 1],
                    alpha=0.6, 
                    c=cmap(idx),
                    edgecolor='black',
                    marker=markers[idx], 
                    label=cl)<jupyter_output><empty_output><jupyter_text>Training logistic regression classifier using the first 2 principal components.<jupyter_code>from sklearn.linear_model import LogisticRegression

lr = LogisticRegression()
lr = lr.fit(X_train_pca, y_train)
plot_decision_regions(X_train_pca, y_train, classifier=lr)
plt.xlabel('PC 1')
plt.ylabel('PC 2')
plt.legend(loc='lower left')
plt.tight_layout()
# plt.savefig('./figures/pca3.png', dpi=300)
plt.show()
plot_decision_regions(X_test_pca, y_test, classifier=lr)
plt.xlabel('PC 1')
plt.ylabel('PC 2')
plt.legend(loc='lower left')
plt.tight_layout()
# plt.savefig('./figures/pca4.png', dpi=300)
plt.show()
lr = LogisticRegression()
lr = lr.fit(X_train_std[:,0:2], y_train)

plot_decision_regions(X_train_std[:,0:2], y_train, classifier=lr)
plt.xlabel('PC 1')
plt.ylabel('PC 2')
plt.legend(loc='lower left')
plt.tight_layout()
# plt.savefig('./figures/pca3.png', dpi=300)
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1487777777777778, 'fraction_numerical': 0.01711111111111111, 'mean_word_length': 4.461771844660194, 'pattern_counts': {'":': 0, '<': 40, '<?xml version=': 0, '>': 43, 'https://': 2, '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': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '232203', 'license': 'MIT', 'n_tokens_mistral': 3418, 'n_tokens_neox': 3174, 'n_words': 751, 'path': 'chap07.ipynb', 'repo_name': 'dev-strender/python-machine-learning-with-scikit-learn', 'size': '9000'}