主成分分析法

目录

  • 主成分分析法

    • 一、主成分分析的理解
    • 二、使用梯度上升法求解PCA
    • 三、求数据的前n个主成分
    • 四、将高维数据向低维数据映射
    • 五、scikit-learn中的PCA
    • 六、对真实数据集MNIST使用PCA
    • 七、使用PCA降噪
    • 八、PCA与人脸识别
    • 我是尾巴：

主成分分析法

主成分分析法：(Principle Component Analysis, PCA)，是一个非监督机器学习算法，主要用于数据降维，通过降维，可以发现便于人们理解的特征，其他应用：可视化和去噪等。

一、主成分分析的理解

? 先假设用数据的两个特征画出散点图，如果我们只保留特征1或者只保留特征2。那么此时就有一个问题，留个哪个特征比较好呢？

? 通过上面对两个特征的映射结果可以发现保留特征1比较好，因为保留特征1，当把所有的点映射到x轴上以后，点和点之间的距离相对较大，也就是说，拥有更高的可区分度，同时还保留着部分映射之前的空间信息。那么如果把点都映射到y轴上，发现点与点距离更近了，这不符合数据原来的空间分布。所以保留特征1相比保留特征2更加合适，但是这是最好的方案吗？

? 也就是说，我们需要找到让这个样本间距最大的轴？那么如何定义样本之间的间距呢？一般我们会使用方差（Variance），Var(x)=\frac{1}{m}\sum_{i=1}^m(x_{i}-\overline{x})^2,找到一个轴，使得样本空间的所有点映射到这个轴的方差最大。

  1. 第一步：将样本进行均值归0，此时$Var(x)=\frac{1}{m}\sum_{i=1}^m(x_{i})^2$，此时计算过程就少一项，这就是为什么要进行样本均值归0。

  1. 第二步：需要定义一个轴的方向w=(w1, w2),使得我们的样本，映射到w以后，有：

Var(X_{project}) = \frac{1}{m}\sum_{i=1}^m(X_{project}^(i) - \overline{x}{project})^2最大。其实括号中的部分是一个向量，更加准确的描述应该是Var(X{project}) = \frac{1}{m}\sum_{i=1}^m\lVertX_{project}^(i) - \overline{x}{project}\lVert^2，因为前面已经去均值，所以，这里只需Var(X{project}) = \frac{1}{m}\sum_{i=1}^m\lVertX_{project}^(i) \lVert^2最大

X^{(i)} \cdot w = \lVert X^{(i)} \lVert \cdot \lVert w \lVert \cdot cos\theta

设w为单位向量，X^{(i)} \cdot w = \lVert X^{(i)} \lVert \cdot cos\theta

最终，X^{(i)} \cdot w = \lVert X_{project}^{(i)} \lVert

所以，最终只需Var(X_{project}) = \frac{1}{m}\sum_{i=1}^m\lVert X^{(i)} \cdot w \lVert^2最大

目标：求w，使得Var(X_{project}) = \frac{1}{m}\sum_{i=1}^m( X^{(i)} \cdot w )^2最大

Var(X_{project}) = \frac{1}{m}\sum_{i=1}^m( X_{1}^{(i)}w_{1} +X_{2}^{(i)}w_{2}+\dots+X_{n}^{(i)}w_{n})^2

Var(X_{project}) = \frac{1}{m}\sum_{i=1}^m( X^{(i)} \cdot w )^2

一个目标函数优化问题，使用梯度上升法解决。

二、使用梯度上升法求解PCA

首先生成一组数据：

import numpy as np
import matplotlib.pyplot as plt

X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0., 10., size=100) 

plt.scatter(X[:, 0], X[:, 1])
plt.show()
def demean(X):
    # 不使用standardscale标准化数据，求解PCA需要方差，只去均值。
    return X - np.mean(X, axis=0) # axis=0表示最终求得的均值是一个行向量，也就是说将每一列求和
x_demean = demean(X)
plt.scatter(x_demean[:,0], x_demean[:,1])
plt.show()
np.mean(x_demean[:,0])
np.mean(x_demean[:,1])
def f(w ,x):
    return np.sum((x.dot(w) ** 2)) / len(x)

def df_math(w, x):
    return x.T.dot(x.dot(w)) * 2. / len(x)

def df_denug(w, x, epsilon=0.0001):
    res = np.empty(len(w))
    for i in range(len(w)):
        w_1 = w.copy()
        w_1[i] += epsilon
        w_2 = w.copy()
        w_2[i] -= epsilon
        res[i] = (f(w_1, x) - f(w_2, x)) / (2 * epsilon)
    return res

def direction(w):
    return w / np.linalg.norm(w)

def gradient_ascent(df, x, init_w, eta, n_iters=1e-4, epsilon=1e-8):

    w = direction(init_w)
    cur_iter = 0
    while cur_iter < n_iters:
        gradient = df(w, x)
        last_w = w
        w = w + eta * gradient
        w = direction(w)
        if (abs(f(w, x) - f(last_w, x)) < epsilon):
            break
        cur_iter += 1
    return w

测试：

init_w = np.random.random(X.shape[1])    # 不能0向量开始
init_w
eta = 0.001
gradient_ascent(df_denug, x_demean, init_w, eta)

输出结果：array([0.84612666, 0.53298187])

w = gradient_ascent(df_math, x_demean, init_w, eta)
plt.scatter(x_demean[:, 0], x_demean[:,1])
plt.plot([0, w[0]*30], [0, w[1]*30], color='r')
plt.show()

这只是一个从2维降到1维的情况，但在实际应用中往往都不这么简单。

三、求数据的前n个主成分

求出第一主成分以后，如何求出下一个主成分？数据进行改变，将数据在第一个主成分上的分量去掉。然后在新的数据上求第一主成分。

第一步：先求出第一主成分

import numpy as np
import matplotlib.pyplot as plt

X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0., 10., size=100) 

def demean(X):
    return X - np.mean(X, axis=0)

X = demean(X)
plt.scatter(X[:,0], X[:,1])
plt.show()

def f(w ,x):
    return np.sum((x.dot(w) ** 2)) / len(x)

def df(w, x):
    return x.T.dot(x.dot(w)) * 2. / len(x)

def direction(w):
    return w / np.linalg.norm(w)

def first_component(df, x, init_w, eta, n_iters=1e-4, epsilon=1e-8):

    w = direction(init_w)
    cur_iter = 0

    while cur_iter < n_iters:
        gradient = df(w, x)
        last_w = w
        w = w + eta * gradient
        w = direction(w)
        if (abs(f(w, x) - f(last_w, x)) < epsilon):
            break
        cur_iter += 1

    return w
init_w = np.random.random(X.shape[1])
init_w
eta = 0.01
w = first_component(df, X, init_w, eta)
w

输出结果: array([0.71849226, 0.69553495])

第二步：去掉第一主成分

X2 = np.empty(X.shape)
for i in range(len(X)):
    X2[i] = X[i] - np.dot(X[i], w) * w
# X2 = X - X.dot(w).reshape(-1, 1) * w

plt.scatter(X2[:,0], X2[:,1])
plt.show()

第三步：求第二主成分

init_w2 = np.random.random(X2.shape[1])
eta = 0.01
w2 = first_component(df, X2, init_w2, eta)
w2

输出结果: array([ 0.98482183, -0.17356834])

第四步：画出主成分

plt.scatter(X2[:,0], X2[:,1])
plt.plot([0, w[0]*30], [0, w[1]*30])
plt.plot([0, w2[0]*30], [0, w2[1]*30])
plt.show()

import numpy as np
import matplotlib.pyplot as plt

X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0., 10., size=100) 

def demean(X):
    return X - np.mean(X, axis=0)

X = demean(X)

def direction(w):
    return w / np.linalg.norm(w)

def df(w, x):
    return x.T.dot(x.dot(w)) * 2. / len(x)

def f(w ,x):
    return np.sum((x.dot(w) ** 2)) / len(x)

ef gradient_ascent(df, x, init_w, eta, n_iters=1e-4, epsilon=1e-8):

    w = direction(init_w)
    cur_iter = 0
    while cur_iter < n_iters:
        gradient = df(w, x)
        last_w = w
        w = w + eta * gradient
        w = direction(w)
        if (abs(f(w, x) - f(last_w, x)) < epsilon):
            break
        cur_iter += 1
    return w

def fisrt_n_component(n, x, eta=0.01, n_iters=1e-4, epsilon=1e-8):
    x_pca = x.copy()
    x_pca = demean(x_pca)
    res = []
    for i in range(n):
        init_w = np.random.random(x_pca.shape[1])
        w = gradient_ascent(df, x_pca, init_w, eta)
        res.append(w)

        x_pca = x_pca - x_pca.dot(w).reshape(-1, 1) * w

    return res    

fisrt_n_component(2, X)

输出结果：[array([0.75978266, 0.65017714]), array([ 0.96416996, -0.26528531])]

四、将高维数据向低维数据映射

? 假设有数据矩阵X，通过主成分分析法得到前k个主成分对应的轴的单位方向向量之后，我们就可以将数据从高维向低维映射，同样也可以从低维在映射回高维空间。

import numpy as np

class PCA(object):

    def __init__(self, n_components):
        assert n_components >= 1, "n_component must be valid"
        self.n_components = n_components
        self.components_ = None

    def fit(self, x, eta=0.01, n_iters=1e-4):
        assert self.n_components <= x.shape[1], "n_components must be greater than the feature number of x"

        def demean(x):
            return x - np.mean(x, axis=0)

        def f(w, x):
            return np.sum((x.dot(w) ** 2)) / len(x)

        def df(w, x):
            return x.T.dot(x.dot(w)) * 2. / len(x)

        def direction(w):
            return w / np.linalg.norm(w)

        def gradient_ascent(df, x, init_w, eta, n_iters=1e-4, epsilon=1e-8):

            w = direction(init_w)
            cur_iter = 0
            while cur_iter < n_iters:
                gradient = df(w, x)
                last_w = w
                w = w + eta * gradient
                w = direction(w)
                if (abs(f(w, x) - f(last_w, x)) < epsilon):
                    break
                cur_iter += 1
            return w

        x_pca = demean(x)
        self.components_ = np.empty(shape=(self.n_components, X.shape[1]))
        for i in range(self.n_components):
            init_w = np.random.random(x_pca.shape[1])
            w = gradient_ascent(df, x_pca, init_w, eta, n_iters)
            self.components_[i, :] = w

            x_pca = x_pca - x_pca.dot(w).reshape(-1, 1) * w

        return self

    def transform(self, x):
        "将给定的x，映射到各个主成分分量中"
        assert x.shape[1] == self.components_.shape[1]
        return x.dot(self.components_.T)

    def inverse_transform(self, x):
        assert x.shape[1] == self.components_.shape[0]
        return x.dot(self.components_)

    def __repr__(self):
        return "PCA(n_components=%d)" % self.n_components
if __name__ == '__main__':
    import numpy as np
    import matplotlib.pyplot as plt

    X = np.empty((100, 2))
    X[:, 0] = np.random.uniform(0., 100., size=100)
    X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0., 10., size=100)

    pca = PCA(n_components=1)
    pca.fit(X)
    print(pca.components_)

    x_reduction = pca.transform(X)
    print(x_reduction.shape)

    x_restore = pca.inverse_transform(x_reduction)
    print(x_restore.shape)

    plt.scatter(X[:, 0], X[:, 1], color='b', alpha=0.5)
    plt.scatter(x_restore[:, 0], x_restore[:, 1], color='r', alpha=0.5)
    plt.show()

? 蓝色的点是随机初始生成，红色的点就是由pca降维之后，在从低维映射到高维的描述。其实完全可以只有一个轴来表示这些红色的点，也就完成了降维。

五、scikit-learn中的PCA

from sklearn.decomposition import PCA
import numpy as np
import matplotlib.pyplot as plt

X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0., 10., size=100) 

def demean(X):
    return X - np.mean(X, axis=0)

X = demean(X)

pca = PCA(n_components=1)
pca.fit(X)
pca.components_

输出结果：array([[-0.7896098, -0.6136093]])

x_reduction = pca.transform(X)
x_reduction.shape
x_restore = pca.inverse_transform(x_reduction)
x_restore.shape
plt.scatter(X[:,0], X[:,1], color='b', alpha=0.5)
plt.scatter(x_restore[:,0], x_restore[:,1], color='r', alpha=0.5)
plt.show()

接下来使用真实的数据集进行测试：

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.model_selection import train_test_split

digits = datasets.load_digits()
x = digits.data
y = digits.target

x_train, x_test, y_train, y_test = train_test_split(x, y, test_size=0.2, random_state=666)
x_train.shape

接下来使用KNN先进行一个简单的测试：

%%time
from sklearn.neighbors import KNeighborsClassifier
knn_clf = KNeighborsClassifier()
knn_clf.fit(x_train, y_train)

输出结果：Wall time: 4.88 ms

knn_clf.score(x_test, y_test)

输出结果：0.9888888888888889

接下来使用降维：

from sklearn.decomposition import PCA
pca = PCA(n_components=2)
pca.fit(x_train)
#PCA(copy=True, iterated_power='auto', n_components=2, random_state=None,svd_solver='auto', tol=0.0, whiten=False)
x_train_reduction = pca.transform(x_train)
x_test_reduction = pca.transform(x_test)
%%time
knn_clf = KNeighborsClassifier()
knn_clf.fit(x_train_reduction, y_train)

输出结果：Wall time: 976 μs 相比降维之前速度快了很多。

knn_clf.score(x_test_reduction, y_test)

输出结果：0.6055555555555555,但是测试结果的准确率大打折扣。这是因为直接从64维降到2维，损失了太多信息。在sklearn中封装了一个函数，可以看出损失的详细信息。

pca.explained_variance_ratio_

输出结果：array([0.1450646 , 0.13714246])，可以看出这两个主成分分别损失了14.51%和13.71%，加和共损失28.22%。显然这不是想要的效果。

那么我们查看一个64个主成分分别损失情况：

pca = PCA(n_components=x_train.shape[1])
pca.fit(x_train)
pca.explained_variance_ratio_

输出结果：

array([1.45064600e-01, 1.37142456e-01, 1.19680004e-01, 8.43768923e-02,
       5.87005941e-02, 5.01797333e-02, 4.34065700e-02, 3.61375740e-02,
       3.39661991e-02, 3.00599249e-02, 2.38906921e-02, 2.29417581e-02,
       1.81335935e-02, 1.78403959e-02, 1.47411385e-02, 1.41290045e-02,
       1.29333094e-02, 1.25283166e-02, 1.01123057e-02, 9.08986879e-03,
       8.98365069e-03, 7.72299807e-03, 7.62541166e-03, 7.09954951e-03,
       6.96433125e-03, 5.84665284e-03, 5.77225779e-03, 5.07732970e-03,
       4.84364707e-03, 4.34595748e-03, 3.73352381e-03, 3.57655938e-03,
       3.30727680e-03, 3.18129431e-03, 3.06969704e-03, 2.89170006e-03,
       2.51205204e-03, 2.27743660e-03, 2.22760483e-03, 2.00065017e-03,
       1.89529684e-03, 1.56877138e-03, 1.42740894e-03, 1.39115781e-03,
       1.20896097e-03, 1.10149976e-03, 9.81702199e-04, 8.82376601e-04,
       5.69898729e-04, 4.10322729e-04, 2.32125043e-04, 8.49807543e-05,
       5.37426557e-05, 5.27990816e-05, 1.03398093e-05, 6.20749843e-06,
       5.03430895e-06, 1.15881302e-06, 1.09689390e-06, 5.51564753e-07,
       5.65215369e-08, 1.78867947e-33, 8.83490979e-34, 8.55393580e-34])

从上面结果可以看出，64个主成分损失从大到小排列，最大14.5%，最小8.55393580e-34可以忽略不计。接下来使用曲线图绘制一下：

plt.plot([i for i in range(x_train.shape[1])],
        [np.sum(pca.explained_variance_ratio_[:i+1]) for i in range(x_train.shape[1])])
plt.show()

从图中，可以大概看出当降到30维的时候，保留信息大概在96%左右，因此可以考虑降到30为左右比较合适，既能保证精度又能保证速度。其实，sklearn中已经封装好了能够实现类似功能。

pca = PCA(0.95)
pca.fit(x_train)
pca.n_components_

输出结果：28,可以发现当损失信息为5%时，可以将数据从64维降到28维。

x_train_reduction = pca.transform(x_train)
x_test_reduction = pca.transform(x_test)
%%time
knn_clf = KNeighborsClassifier()
knn_clf.fit(x_train_reduction, y_train)

输出结果：Wall time: 2.93 ms

knn_clf.score(x_test_reduction, y_test)

输出结果：0.9833333333333333

对比降维前后：

  1. 时间上，降维前：Wall time: 4.88 ms，降维后：Wall time: 2.93 ms
  2. 精度上，降维前：0.9888888888888889，降维后0.9833333333333333

通过以上对比，可以发现在精度损失0.5%的情况下速度能提高一倍。这是一个比较的维度，如果维度很大，效果会更加优秀。

在整个测试的过程中，曾尝试把数据降到2维，精度损失很大，但并不能说明这样做没有意义，通常情况下，我们会将高维数据降到3维或者2维进行可视化。

pca = PCA(n_components=2)
pca.fit(x)
x_reduction = pca.transform(x)
x_reduction.shape
for i in range(10):
    plt.scatter(x_reduction[y==i, 0], x_reduction[y==i, 1], alpha=0.8)
plt.show()

通过对数据进行可视化，可以发现有蓝色、红色、紫色这些点在降到2维的情况下，依旧可以实现很高的识别度。这说明这些样本在原来的高维空间中差异就比较大。所以，一般会将数据进行可视化，进行一些简单的分析。

六、对真实数据集MNIST使用PCA

首先下载MNIST数据集的时候有个坑？

import numpy as np
from sklearn.datasets import fetch_mldata

mnist = fetch_mldata("MNIST original", data_home='./datasets')

如果直接运行会一直报错，反正就是下载不下来，但是会在datasets目录下生成一个文件夹，./datasets/mldata。所以需要下载数据集放在这个文件夹下，再次运行就可以了。

数据集地址：

链接：https://pan.baidu.com/s/15k6Sykf5TICN40KeEc6xgQ
提取码：ujrq

mnist
{'DESCR': 'mldata.org dataset: mnist-original',
 'COL_NAMES': ['label', 'data'],
 'target': array([0., 0., 0., ..., 9., 9., 9.]),
 'data': array([[0, 0, 0, ..., 0, 0, 0],
        [0, 0, 0, ..., 0, 0, 0],
        [0, 0, 0, ..., 0, 0, 0],
        ...,
        [0, 0, 0, ..., 0, 0, 0],
        [0, 0, 0, ..., 0, 0, 0],
        [0, 0, 0, ..., 0, 0, 0]], dtype=uint8)}
x, y = mnist['data'], mnist['target']
x.shape
y.shape

输出结果：(70000, 784) (70000,)

x_train = np.array(x[:60000], dtype=float)
y_train = np.array(y[:60000], dtype=float)
x_test = np.array(x[60000:], dtype=float)
x_test = np.array(y[60000:], dtype=float)
from sklearn.neighbors import KNeighborsClassifier

knn_clf = KNeighborsClassifier()
%time knn_clf.fit(x_train, y_train)

输出结果：

Wall time: 28 s
KNeighborsClassifier(algorithm='auto', leaf_size=30, metric='minkowski',
           metric_params=None, n_jobs=None, n_neighbors=5, p=2,
           weights='uniform')
%time knn_clf.score(x_test, y_test)

输出结果：

Wall time: 11min 6s
0.9688

在这里，值得一提的是，之前讨论knn的时候，由于相比较样本之间距离，所以通常会对数据进行归一化。但是这里没用使用standardscale，因为这是图像。因为图像的像素点都在0-255，所以进行距离的比较是有意的，只有当尺度不一致时，才会考虑使用归一化。

接下来使用PCA进行降维：

from sklearn.decomposition import PCA

pca = PCA(0.9)
pca.fit(x_train)
x_train_reduction = pca.transform(x_train)
x_train.shape

输出结果：(60000, 87) 可以发现通过降维实现了只有87维代替原来的784维，能够保留90的信息。

from sklearn.neighbors import KNeighborsClassifier

knn_clf = KNeighborsClassifier()
%time knn_clf.fit(x_train_reduction, y_train)

输出结果：Wall time: 436 ms 降维前使用了28s完成现在只需436ms。

x_test_reduction = pca.transform(x_test)
%time knn_clf.score(x_test_reduction, y_test)
Wall time: 1min 17s
0.9728

将降维前后进行对比：

  1. 时间上：

    • 训练时间-降维前：28 s，降维后：436 ms
    • 测试时间-降维前：11min 6s，降维后：1min 17s
  2. 精度上：降维前：0.9688，降维后：0.9728

经过上面的对比，发现经过降维之后速度提高了很多很多，而且准确率反而高了，这是为什么呢？这是因为其实降维还有另外一个作用，就是去噪。

七、使用PCA降噪

import numpy as np
import matplotlib.pyplot as plt

X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 5, size=100)
plt.scatter(X[:,0], X[:,1])
plt.show()
from sklearn.decomposition import PCA
pca = PCA(n_components=1)
pca.fit(X)
x_reduction = pca.transform(X)
x_restore = pca.inverse_transform(x_reduction)
plt.scatter(x_restore[:, 0], x_restore[:, 1])
plt.show()

通过这个例子，发现pca还可以降噪，将高维数据降低到低维数据，这个过程中也可以理解成降低了维度，丢失了信息，同时也去除了噪音。为了更加直观，接下来使用手写数字识别数据集。

from sklearn import datasets

digits = datasets.load_digits()
x = digits.data
y = digits.target
# 手动添加噪声
noisy_digits = x + np.random.normal(0, 4, size=x.shape)
# 取出部分数据做可视化
example_digits = noisy_digits[y==0,:][:10]
for num in range(1, 10):
    x_num = noisy_digits[y==num,:][:10]
    example_digits = np.vstack([example_digits, x_num])
example_digits.shape
def plot_digits(data):
    fig, axes = plt.subplots(10, 10, figsize=(10,10),
                            subplot_kw={'xticks':[], 'yticks':[]},
                            gridspec_kw=dict(hspace=0.1, wspace=0.1))
    for i, ax in enumerate(axes.flat):
        ax.imshow(data[i].reshape(8, 8),
                 cmap='binary', interpolation='nearest',
                 clim=(0, 16))
    plt.show()

 plot_digits(example_digits)

接下来使用pca去噪，

from sklearn.decomposition import PCA

pca = PCA(0.5)
pca.fit(noisy_digits)
pca.n_components_

输出结果：12，可以看出保留50%的信息，只需要12维。

components = pca.transform(example_digits)
filterd_digits = pca.inverse_transform(components)
plot_digits(filterd_digits)

通过这个例子，可以对pca降维有一个感性的认识。

八、PCA与人脸识别

假设X一张人脸图像，W是钱k个主成分组成的矩阵，其实可以把W看做X的一些比较明显的特征。这就是特征脸（Eigenface）。

首先下载数据集，sklearn中的lfw_people数据集：

import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_lfw_people

faces = fetch_lfw_people()

这里依旧有个小坑，就是数据下载不下来，不过没事，我们可以手动下载，复制这个网址 https://ndownloader.figshare.com/files/5976015 到迅雷中，下载到完整的lfwfunneded.tgz文件。并把这个文件放在C:\Users\Mogul\scikit_learn_data\lfw_home路径下，解压。然后再次运行就可以了。

faces
faces.keys()
faces.data.shape
#(13233, 2914)，数据集中一共13233个样本，每个样本维度2914
faces.images.shape
# (13233, 62, 47)，每个图片的大小（62， 47）
random_indexes = np.random.permutation(len(faces.data))
# 将数据集的顺序打乱
X = faces.data[random_indexes]
# 从中取除前36个
example_faces = X[:36, :]
example_faces.shape
def plot_faces(faces):
    fig, axes = plt.subplots(6, 6, figsize=(10, 10),
                            subplot_kw={'xticks':[], 'yticks':[]},
                            gridspec_kw=dict(hspace=0.1, wspace=0.1))
    for i, ax in enumerate(axes.flat):
        ax.imshow(faces[i].reshape(62, 47), cmap='bone')

    plt.show()
plot_faces(example_faces)

faces.target_names
# array(['AJ Cook', 'AJ Lamas', 'Aaron Eckhart', ..., 'Zumrati Juma',
       'Zurab Tsereteli', 'Zydrunas Ilgauskas'], dtype='<U35')
len(faces.target_names)
# 5749,数据集中一共有5749个人物，一共有13233张图片，差不多每个人2张。
# 但实际情况并不是这样，因为样本并不均衡，有些可能只有一张，有些可能很多张。
%%time
from sklearn.decomposition import PCA
pca = PCA(svd_solver='randomized')
pca.fit(X)

输出结果：Wall time: 18.4 s

pca.components_.shape
# (2914, 2914)
# 取出前36个特征脸进行可视化
plot_faces(pca.components_[:36, :])

刚刚提到这个数据集的样本分布是不均衡的，那么接下来看看具体：

faces2 = fetch_lfw_people(min_faces_per_person=60)
faces2.data.shape
# (1348, 2914)
faces2.target_names

输出结果：

array(['Ariel Sharon', 'Colin Powell', 'Donald Rumsfeld', 'George W Bush',
       'Gerhard Schroeder', 'Hugo Chavez', 'Junichiro Koizumi',
       'Tony Blair'], dtype='<U17')

根据输出结果可以发现，最少有60张图像的数据集一共有1348张图片，只有8个人的图像大于60，使用这样的数据进行人脸识别相对就比较合理了。

我是尾巴：

其实在PCA的背后是强大的数学支撑。在这里使用的梯度上升法对PCA进行求解，其实它是可以通过数学推导出相应的公式，具体的数学解释，以后会在写一篇。

本次推荐：

轻松玩转PDF

你是如何强迫自己不断学习提升的？

看过更大的世界，更优秀的人，就再也不甘心留在原地，不甘心就是动力！

继续加油~

原文地址：https://www.cnblogs.com/zhangkanghui/p/11331276.html

时间： 2024-10-11 06:14:44

主成分分析法的相关文章

主成分分析法原理及其python实现

主成分分析法原理及其python实现 前言: 这片文章主要参考了Andrew Ng的Machine Learning课程讲义,我进行了翻译,并配上了一个python演示demo加深理解. 本文主要介绍一种降维算法,主成分分析法,Principal Components Analysis,简称PCA,这种方法的目标是找到一个数据近似集中的子空间,至于如何找到这个子空间,下文会给出详细的介绍,PCA比其他降维算法更加直接,只需要进行一次特征向量的计算即可.(在Matlab,python,R中这个可以

[数学模型]主成分分析法python实现

def pca(dataMat, topNfeat=9999999): #数据矩阵, 输出前topNfeat个feat meanVals = mean(dataMat, axis=0) # 计算平均值 meanRemoved = dataMat - meanVals covMat = cov(meanRemoved, rowvar=0) #计算协方差矩阵 eigVals,eigVects = linalg.eig(mat(covMat)) #特征值, eigValInd = argsort(ei

算法理解之主成分分析法

主成分分析法就是通过正交变换将存在相关性的原始变量变量转换成不相关的新变量,将其中贡献度低的变量舍弃掉,贡献度可以理解为变量的方差,方差越大,贡献度越高,正交变换前后变量存在的信息量是相同的,只是把更多的信息集中起来,舍弃存在少量信息的变量,达到降维的目的. 运用主成分的前提是变量之间存在相关性 主成分分析法的实质的正交变换,正交变换实质就是矩阵的对角化,使得协方差矩阵为对角矩阵,各个变量之间的相关性为0. 如何确定主成分数量,累计贡献度达到90%左右 如何确定主成分所代表的的含义:通过找到荷载

降维之主成分分析法（PCA）

一.主成分分析法的思想 我们在研究某些问题时,需要处理带有很多变量的数据,比如研究房价的影响因素,需要考虑的变量有物价水平.土地价格.利率.就业率.城市化率等.变量和数据很多,但是可能存在噪音和冗余,因为这些变量中有些是相关的,那么就可以从相关的变量中选择一个,或者将几个变量综合为一个变量,作为代表.用少数变量来代表所有的变量,用来解释所要研究的问题,就能从化繁为简,抓住关键,这也就是降维的思想. 主成分分析法(Principal Component Analysis,PCA)就是一种运用线性代

PCA(principal component analysis)主成分分析法

<Aggregating local descriptors into a compact image representation>论文笔记 在论文中,提取到VLAD特征后,要对特征向量进行PCA降维,就是用一个大小为D' * D的矩阵M,对VLAD特征向量x做变换,降维后的vector是x' = Mx,x'的大小是D'维.矩阵M是由原样本的协方差矩阵的D'个特征向量构成. 为什么M要是特征向量的矩阵呢? 根据PRML中的内容,理解如下: 1,Maxinum Variance Formula

R语言学习之主成分分析法的R实践

主成分分析R软件实现程序(一): >d=read.table("clipboard",header=T) #从剪贴板读取数据 >sd=scale(d)  #对数据进行标准化处理 >sd  #输出标准化后的数据和属性信息,把标准化的数据拷贝到剪贴板备用 >d=read.table("clipboard",header=T)  #从剪贴板读取标准化数据 >pca=princomp(d,cor=T)  #主成分分析函数 >screepl

【cs229-Lecture14】主成分分析法

本节课内容: 因子分析 ---因子分析中的EM步骤的推导过程 主成份分析:有效地降低维度的方法 因子分析 混合高斯模型的问题 接下来讨论因子分析模型,在介绍因子分析模型之前,先看高斯分布的另一种写法,该写法是推导因子分析模型的基础. 高斯分布的矩阵写法 因子分析模型 因子分析模型的推导 EM 求解参数 PCA(Principal  Components  Analysis,  主成分分析),也是一种降维方法 主要介绍 PCA(Principal  Components  Analysis,  主

主成分分析法(PCA)答疑

问:为什么要去均值? 1.我认为归一化的表述并不太准确,按统计的一般说法,叫标准化.数据的标准化过程是减去均值并除以标准差.而归一化仅包含除以标准差的意思或者类似做法.2.做标准化的原因是:减去均值等同于坐标的移动,把原始数据点的重心移到和原点重合,这样利于很多表达,比如数据的协方差矩阵可以写成XX',若没有减去均值,则XX‘后面还要减去一些东西(还不明白可以参考多元统计分析的书).除以标准差是为了统一并消除量纲.一个矩阵中有多个向量,有些可能表示了长度,有些表示了重量,除以标准差,才能让它们仅

主成分分析法PCA原理

PCA(Principal Component Analysis)是一种常用的数据分析方法.PCA通过线性变换将原始数据变换为一组各维度线性无关的表示,可用于提取数据的主要特征分量,常用于高维数据的降维.网上关于PCA的文章有很多,但是大多数只描述了PCA的分析过程,而没有讲述其中的原理.这篇文章的目的是介绍PCA的基本数学原理,帮助读者了解PCA的工作机制是什么. 当然我并不打算把文章写成纯数学文章,而是希望用直观和易懂的方式叙述PCA的数学原理,所以整个文章不会引入严格的数学推导.希望读者在