These figures aid in illustrating how a point cloud can be very flat in one direction?which is where PCA comes in to choose a direction that is not flat.
1 2 3 4 5 6 7 8 9 10 11 12 13 | print (__doc__) # Authors: Gael Varoquaux # Jaques Grobler # Kevin Hughes # License: BSD 3 clause from sklearn.decomposition import PCA from mpl_toolkits.mplot3d import Axes3D import numpy as np import matplotlib.pyplot as plt from scipy import stats |
Create the data
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 | e = np.exp( 1 ) np.random.seed( 4 ) def pdf(x): return 0.5 * (stats.norm(scale = 0.25 / e).pdf(x) + stats.norm(scale = 4 / e).pdf(x)) y = np.random.normal(scale = 0.5 , size = ( 30000 )) x = np.random.normal(scale = 0.5 , size = ( 30000 )) z = np.random.normal(scale = 0.1 , size = len (x)) density = pdf(x) * pdf(y) pdf_z = pdf( 5 * z) density * = pdf_z a = x + y b = 2 * y c = a - b + z norm = np.sqrt(a.var() + b.var()) a / = norm b / = norm |
Plot the figures
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 | def plot_figs(fig_num, elev, azim): fig = plt.figure(fig_num, figsize = ( 4 , 3 )) plt.clf() ax = Axes3D(fig, rect = [ 0 , 0 , . 95 , 1 ], elev = elev, azim = azim) ax.scatter(a[:: 10 ], b[:: 10 ], c[:: 10 ], c = density[:: 10 ], marker = '+' , alpha = . 4 ) Y = np.c_[a, b, c] # Using SciPy's SVD, this would be: # _, pca_score, V = scipy.linalg.svd(Y, full_matrices=False) pca = PCA(n_components = 3 ) pca.fit(Y) pca_score = pca.explained_variance_ratio_ V = pca.components_ x_pca_axis, y_pca_axis, z_pca_axis = V.T * pca_score / pca_score. min () x_pca_axis, y_pca_axis, z_pca_axis = 3 * V.T x_pca_plane = np.r_[x_pca_axis[: 2 ], - x_pca_axis[ 1 :: - 1 ]] y_pca_plane = np.r_[y_pca_axis[: 2 ], - y_pca_axis[ 1 :: - 1 ]] z_pca_plane = np.r_[z_pca_axis[: 2 ], - z_pca_axis[ 1 :: - 1 ]] x_pca_plane.shape = ( 2 , 2 ) y_pca_plane.shape = ( 2 , 2 ) z_pca_plane.shape = ( 2 , 2 ) ax.plot_surface(x_pca_plane, y_pca_plane, z_pca_plane) ax.w_xaxis.set_ticklabels([]) ax.w_yaxis.set_ticklabels([]) ax.w_zaxis.set_ticklabels([]) elev = - 40 azim = - 80 plot_figs( 1 , elev, azim) elev = 30 azim = 20 plot_figs( 2 , elev, azim) plt.show() |
Total running time of the script: (0 minutes 0.278 seconds)
Download Python source code:
plot_pca_3d.py
Download IPython notebook:
plot_pca_3d.ipynb
Please login to continue.