分解成分信号

一、主成分分析

1.1 主成分分析

在多元统计分析中,主成分分析(Principal components analysis,PCA)是一种分析、简化数据集的技术。主成分分析经常用于减少数据集的维数,同时保持数据集中的对方差贡献最大的特征。这是通过保留低阶主成分,忽略高阶主成分做到的。这样低阶成分往往能够保留住数据的最重要方面。但是,这也不是一定的,要视具体应用而定。由于主成分分析依赖所给数据,所以数据的准确性对分析结果影响很大。

主成分分析由卡尔·皮尔逊于1901年发明,用于分析数据及建立数理模型。其方法主要是通过对协方差矩阵进行特征分解,以得出数据的主成分(即特征向量)与它们的权值(即特征值)。PCA是最简单的以特征量分析多元统计分布的方法。其结果可以理解为对原数据中的方差做出解释:哪一个方向上的数据值对方差的影响最大?换而言之,PCA提供了一种降低数据维度的有效办法;如果分析者在原数据中除掉最小的特征值所对应的成分,那么所得的低维度数据必定是最优化的(也即,这样降低维度必定是失去讯息最少的方法)。主成分分析在分析复杂数据时尤为有用,比如人脸识别。

PCA是最简单的以特征量分析多元统计分布的方法。通常情况下,这种运算可以被看作是揭露数据的内部结构,从而更好的解释数据的变量的方法。如果一个多元数据集能够在一个高维数据空间坐标系中被显现出来,那么PCA就能够提供一幅比较低维度的图像,这幅图像即为在讯息最多的点上原对象的一个‘投影’。这样就可以利用少量的主成分使得数据的维度降低了。

PCA跟因子分析密切相关,并且已经有很多混合这两种分析的统计包。而真实要素分析则是假定底层结构,求得微小差异矩阵的特征向量。

主成分分析(Principal Component Analysis)是最常用的一种降维方法,是利用在正交属性空间上的样本点,用一个超平面对所有样本进行恰当的表达。

如果这里我们可以用直线的斜率和截距(k,b)来表示超平面的一些点,点在直线上的投影点(这时可以把直线看作一维的)来表示。若存在这样的超平面,则应具有的性质为

假设D=\{x_1,x_2,\dots,x_m\}为样本集,大小为m.

最近重构性: 假定对数据样本进行了中心化,即\Sigma_ix_i=0,投影变换后得到的新坐标系为\{\omega_1,\omega_2,\dots,\omega_d\},其中\omega_i是标准正交基向量:\|\omega_i\|_2=1,\omega_i^T\omega_j=0(i\neq j),若丢弃新坐标系中的部分向量,则可将维度降为d'<d,假设现在丢弃后的坐标系为W=\{\omega_1,\omega_2,\dots,\omega_{d'}\},则样本点在低维坐标系中的投影为z_i=\{z_{i1};z_{i2};\dots;z_{id'}\},其中z_{ij}=\omega_j^Tx_ix_i在新坐标系下第j维的坐标,若基于z_i来重构x_i,则得到\hat{x}_i=\Sigma_{j=1}^{d'}z_{ij}\omega_j

对整个训练集,原样本点x_i与基于投影重构的样本点\hat{x}_i之间的距离为

\begin{aligned} \sum\limits_{i=1}^{m}||\sum\limits_{j=1}^{d'}z_{ij}\omega_j-x_i||_2^2 & = \sum\limits_{i=1}^{m}||Wz_i-x_i||_2^2 \\ & = \sum\limits_{i=1}^mz_i^Tz_i-2\sum\limits_{i=1}^mz_i^TW^Tx_i+const(\sum x_i^Tx_i) \\ & = -tr(W^T(\sum\limits_{i=1}^mx_ix_i^T)W) \end{aligned}

最后一步利用了tr(AB)=tr(BA)z_i=W^Tx_i。根据最近重构性,可得优化目标函数

\begin{aligned} &\min\limits_{W}\quad -tr(W^TXX^TW) \\ &s.t.\quad W^TW=I \end{aligned}

其中\omega_j是标准正交基,\sum_ix_ix_i^T是协方差矩阵。

最大可分性: 样本点x_i在新空间上的投影为W^Tx_i,若所有样本点能近可能分开,则应该使投影后样本的方差最大,即\sum_ix_i^TWW^Tx_i,由x_i^TWW^Tx_i=tr(W^Tx_ix_i^TW),可知优化目标函数为

\begin{aligned} &\max\limits_W \quad tr(W^TXX^TW) \\ &s.t. \quad W^TW=I \end{aligned}

求解优化函数:

\begin{aligned} &\max\limits_W \quad tr(W^TXX^TW) \\ &s.t. \quad W^TW=I \end{aligned}

假设d'=1,构造Lagrange函数,可以得到

f(\omega)=\omega^TXX^T\omega+\lambda(\omega^T\omega)

\omega求导,得到XX^T\omega=\lambda\omega,可以知道\omegaXX^T的特征向量。

假设d'>1,构造Lagrange函数,可以得到

f(\omega_1,\cdots,\omega_{d'})=\sum\limits_{i=1}^{d'}\omega_i^TXX^T\omega_i+\lambda_{ij}(\omega_i^T\omega_j-\delta_{ij})

\omega_i求导,得到

XX^T\omega_i=\lambda_{ii}\omega_i+(\lambda_{ij}+\lambda_{ji})\omega_j

从PCA的目标,我们知道,我们要确定一个超平面来表示x_1,\cdots,x_n,但是表示这个超平面的基向量可以不唯一,但是这个超平面必须是唯一的,由上面公式可知表示这个超平面的基向量可以互相转换,而且XX^T在这个超平面内恰有d'个特征向量,所以我们可以用XX^T的特征向量来表示这个超平面。即优化函数的解为

XX^T\omega=\lambda \omega

python实现:

输入:样本集D=\{x_1,\dots ,x_m\},低维空间维数d'

算法过程:

  1. 对所有样本进行中心化x_i\leftarrow \frac{1}{m}\sum_{i=1}^mx_i;
  2. 计算样本的协方差XX^T;
  3. 对协方差XX^T做特征值分解;
  4. 取最大的d'个特征值所对应的特征向量\omega_1,\omega_2,\dots,\omega_{d'};

输出:投影矩阵W=\{\omega_1,\dots,\omega_{d'}\}

import numpy as np
import matplotlib.pyplot as plt

# PCA获取主特征向量和均值
def PCA(date,args=1):
    mean=np.sum(date,1)/mat.shape[1];
    temp_mat=date-np.outer(mean,np.ones(mat.shape[1]));    # 中心化的矩阵
    date_shape=temp_mat.shape
    if date_shape[0]>date_shape[1]:
        covar_mat=np.dot(temp_mat.T,temp_mat);      # 计算协方差
        eigval,eigvec = np.linalg.eig(covar_mat);   # 计算特征值和特征向量
        eigvec=np.dot(covar_mat,eigvec)
    else:
        covar_mat=np.dot(temp_mat,temp_mat.T)
        eigval,eigvec = np.linalg.eig(covar_mat);

    indices=np.argsort(eigval)[::-1];           # 得到特征值降序排列的下标
    eigval=eigval[indices];                     # 得到排序的特征值
    eigvec=eigvec[:,indices];                     # 得到排序的特征向量
    if args>=1:
        return eigvec[:,0:int(args)],mean
    else:
        temp_sum,d,eig_sum=0,0,np.sum(eigval)
        while temp_sum/eig_sum<args:
            temp_sum=temp_sum+eigval[d]
            d=d+1
        return eigvec[:,0:d],mean;

# 降维获取相应的坐标
def dime_reduction(date,eigvec,mean):
    return np.dot(vec.T,date-np.outer(mean,np.ones(date.shape[1])))

# 测试数据
x,y=np.random.multivariate_normal([10,10],[[2,1],[1,2]],1000).T
mat=np.array([list(x),list(y)]);
vec,mean=PCA(mat,1);
indice=dime_reduction(mat,vec,mean)
PCA_date=np.outer(vec,indice)+np.outer(mean,np.ones(mat.shape[1]))
plt.plot(mat[0],mat[1],'r.')
plt.plot(PCA_date[0],PCA_date[1],'g')
plt.show()

其中d'还可以由其他方法获得,如设置重构阈值t=95%,s使得\frac{\sum_{i=1}^{d'}\lambda_i}{\sum_{i=1}^d\lambda_i}>t成立的最小d'.

注意: 降维导致了舍弃了一部分信息,这往往是必要的,一是舍弃这部分信息后,使得样本的采样密度增大(降维的重要动机);二是最小的特征值所对应的特征值往往与噪声有关,将他们舍弃在一定程度上起到去噪的作用.

PCA的数学定义是:一个正交化线性变换,把数据变换到一个新的坐标系统中,使得这一数据的任何投影的第一大方差在第一个坐标(称为第一主成分)上,第二大方差在第二个坐标(第二主成分)上,依次类推。

定义一个n×m的矩阵,X^T为去平均值(以平均值为中心移动至原点)的数据,其行为数据样本,列为数据类别(注意:这里定义的是X^T,而不是X)。则X的奇异值分解为X = W\Sigma V^T,其中m×m矩阵WXX^T的本征矢量矩阵,\Sigmam×n的非负矩形对角矩阵,Vn×nX^TX的本征矢量矩阵。据此,

\begin{aligned} {\boldsymbol {Y}}^{\top }&={\boldsymbol {X}}^{\top }{\boldsymbol {W}}\\ &={\boldsymbol {V}}{\boldsymbol {\Sigma }}^{\top }{\boldsymbol {W}}^{\top }{\boldsymbol {W}}\\ &={\boldsymbol {V}}{\boldsymbol {\Sigma }}^{\top }\end{aligned}

m<n-1时,V在通常情况下不是唯一定义的,而Y则是唯一定义的。W是一个正交矩阵,Y^TW^T=X^T,且Y^T的第一列由第一主成分组成,第二列由第二主成分组成,依此类推。

为了得到一种降低数据维度的有效办法,我们可以利用W_LX映射到一个只应用前面L个向量的低维空间中去:

\mathbf {Y} =\mathbf {W_{L}} ^{\top }\mathbf {X} =\mathbf {\Sigma _{L}} \mathbf {V} ^{\top }

其中\mathbf {\Sigma _{L}} =\mathbf {I} _{L\times m}\mathbf {\Sigma },且 \mathbf {I} _{L\times m}L\times m的单位矩阵。

X的单向量矩阵W相当于协方差矩阵的本征矢量C = X X^T, \mathbf {X} \mathbf {X} ^{\top }=\mathbf {W} \mathbf {\Sigma } \mathbf {\Sigma } ^{\top }\mathbf {W} ^{\top }

在欧几里得空间给定一组点数,第一主成分对应于通过多维空间平均点的一条线,同时保证各个点到这条直线距离的平方和最小。去除掉第一主成分后,用同样的方法得到第二主成分。依此类推。在\Sigma中的奇异值均为矩阵XX^T的本征值的平方根。每一个本征值都与跟它们相关的方差是成正比的,而且所有本征值的总和等于所有点到它们的多维空间平均点距离的平方和。PCA提供了一种降低维度的有效办法,本质上,它利用正交变换将围绕平均点的点集中尽可能多的变量投影到第一维中去,因此,降低维度必定是失去讯息最少的方法。PCA具有保持子空间拥有最大方差的最优正交变换的特性。然而,当与离散余弦变换相比时,它需要更大的计算需求代价。非线性降维技术相对于PCA来说则需要更高的计算要求。

PCA对变量的缩放很敏感。如果我们只有两个变量,而且它们具有相同的样本方差,并且成正相关,那么PCA将涉及两个变量的主成分的旋转。但是,如果把第一个变量的所有值都乘以100,那么第一主成分就几乎和这个变量一样,另一个变量只提供了很小的贡献,第二主成分也将和第二个原始变量几乎一致。这就意味着当不同的变量代表不同的单位(如温度和质量)时,PCA是一种比较武断的分析方法。但是在Pearson的题为 "On Lines and Planes of Closest Fit to Systems of Points in Space"的原始文件里,是假设在欧几里得空间里不考虑这些。一种使PCA不那么武断的方法是使用变量缩放以得到单位方差。

通常,为了确保第一主成分描述的是最大方差的方向,我们会使用平均减法进行主成分分析。如果不执行平均减法,第一主成分有可能或多或少的对应于数据的平均值。另外,为了找到近似数据的最小均方误差,我们必须选取一个零均值。

假设零经验均值,数据集X 的主成分w_1可以被定义为:

\mathbf{w}_1=\underset{\Vert \mathbf{w} \Vert = 1}{\operatorname{\arg\,max}}\,\operatorname{Var}\{ \mathbf{w}^\top \mathbf{X} \} = \underset{\Vert \mathbf{w} \Vert = 1}{\operatorname{\arg\,max}}\,E\left\{ \left( \mathbf{w}^\top \mathbf{X}\right)^2 \right\}

为了得到第k个主成分,必须先从X中减去前面的k-1个主成分:

\mathbf{\hat{X}}_{k-1}=\mathbf{X}-\sum_{i = 1}^{k - 1}\mathbf{w}_i \mathbf{w}_i^\top \mathbf{X}

然后把求得的第k个主成分带入数据集,得到新的数据集,继续寻找主成分。

\mathbf{w}_k= \underset{\Vert \mathbf{w} \Vert = 1}{\operatorname{arg\,max}}\,E\left\{\left( \mathbf{w}^\top \mathbf{\hat{X}}_{k-1}\right)^2 \right\}.

PCA相当于在气象学中使用的经验正交函数(EOF),同时也类似于一个线性隐层神经网络。隐含层K个神经元的权重向量收敛后,将形成一个由前K个主成分跨越空间的基础。但是与PCA不同的是,这种技术并不一定会产生正交向量。

PCA是一种很流行且主要的的模式识别技术。然而,它并不能最优化类别可分离性 。另一种不考虑这一点的方法是线性判别分析。

class sklearn.decomposition.PCA(n_components=None, *, copy=True, whiten=False, svd_solver='auto')
    """
    Parameters
        n_componentsint, float, None or str
            Number of components to keep. if n_components is not set all components are kept:
                n_components == min(n_samples, n_features)
                If n_components == 'mle' and svd_solver == 'full', Minka’s MLE is used to guess the dimension. Use of n_components == 'mle' will interpret svd_solver == 'auto' as svd_solver == 'full'.
                If 0 < n_components < 1 and svd_solver == 'full', select the number of components such that the amount of variance that needs to be explained is greater than the percentage specified by n_components.
                If svd_solver == 'arpack', the number of components must be strictly less than the minimum of n_features and n_samples.
        copy bool, default=True
            If False, data passed to fit are overwritten and running fit(X).transform(X) will not yield the expected results, use fit_transform(X) instead.
        whiten bool, optional (default False)
            When True (False by default) the components_ vectors are multiplied by the square root of n_samples and then divided by the singular values to ensure uncorrelated outputs with unit component-wise variances.
            Whitening will remove some information from the transformed signal (the relative variance scales of the components) but can sometime improve the predictive accuracy of the downstream estimators by making their data respect some hard-wired assumptions.
        svd_solver str {‘auto’, ‘full’, ‘arpack’, ‘randomized’}
            If auto : The solver is selected by a default policy based on X.shape and n_components: if the input data is larger than 500x500 and the number of components to extract is lower than 80% of the smallest dimension of the data, then the more efficient ‘randomized’ method is enabled. Otherwise the exact full SVD is computed and optionally truncated afterwards.
            If full :run exact full SVD calling the standard LAPACK solver via scipy.linalg.svd and select the components by postprocessing
            If arpack : run SVD truncated to n_components calling ARPACK solver via scipy.sparse.linalg.svds. It requires strictly 0 < n_components < min(X.shape)
            If randomized : run randomized SVD by the method of Halko et al.
    Attributes
        components_array, shape (n_components, n_features)
            Principal axes in feature space, representing the directions of maximum variance in the data. The components are sorted by explained_variance_.
        explained_variance_array, shape (n_components,)
            The amount of variance explained by each of the selected components.
            Equal to n_components largest eigenvalues of the covariance matrix of X.
    Examples
        >>> import numpy as np
        >>> from sklearn.decomposition import PCA
        >>> X = np.array([[-1, -1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
        >>> pca = PCA(n_components=2)
        >>> pca.fit(X)
        PCA(n_components=2)
        >>> print(pca.explained_variance_ratio_)
        [0.9924... 0.0075...]
        >>> print(pca.singular_values_)
        [6.30061... 0.54980...]
        >>>
        >>> pca = PCA(n_components=2, svd_solver='full')
        >>> pca.fit(X)
        PCA(n_components=2, svd_solver='full')
        >>> print(pca.explained_variance_ratio_)
        [0.9924... 0.00755...]
        >>> print(pca.singular_values_)
        [6.30061... 0.54980...]
        >>>
        >>> pca = PCA(n_components=1, svd_solver='arpack')
        >>> pca.fit(X)
        PCA(n_components=1, svd_solver='arpack')
        >>> print(pca.explained_variance_ratio_)
        [0.99244...]
        >>> print(pca.singular_values_)
        [6.30061...]
    """

方法

方法 说明
fit(self, X[, y]) Fit the model with X.
fit_transform(self, X[, y]) Fit the model with X and apply the dimensionality reduction on X.
get_covariance(self) Compute data covariance with the generative model.
get_params(self[, deep]) Get parameters for this estimator.
get_precision(self) Compute data precision matrix with the generative model.
inverse_transform(self, X) Transform data back to its original space.
score(self, X[, y]) Return the average log-likelihood of all samples.
score_samples(self, X) Return the log-likelihood of each sample.
set_params(self, **params) Set the parameters of this estimator.
transform(self, X) Apply dimensionality reduction to X.

1.2 增量PCA(Incremental PCA)

PCA对象非常有用,但针对大型数据集的应用,仍然具有一定的限制。最大的限制是PCA仅支持批处理,这意味着所有要处理的数据必须放在主内存。IncrementalPCA对象使用不同的处理形式,即允许部分计算以小型批处理方式处理数据的方法进行,而得到和PCA算法差不多的结果。 IncrementalPCA可以通过以下方式实现核外(out-of-core)主成分分析:

基于从本地硬盘或网络数据库中连续获取的数据块之上,使用partial_fit方法

在 memory mapped file(通过numpy.memmap创建)上使用fit方法。

IncrementalPCA 类为了增量式的更新explainedvariance_ratio,仅需要存储估计出的分量和噪声方差。 在应用SVD之前,IncrementalPCA就像PCA一样,为每个特征聚集而不是缩放输入数据。

class sklearn.decomposition.IncrementalPCA(n_components=None, *, whiten=False, copy=True, batch_size=None)
    """
    Incremental principal components analysis (IPCA).

    Linear dimensionality reduction using Singular Value Decomposition of the data, keeping only the most significant singular vectors to project the data to a lower dimensional space. The input data is centered but not scaled for each feature before applying the SVD.

    Depending on the size of the input data, this algorithm can be much more memory efficient than a PCA, and allows sparse input.

    Parameters
        n_componentsint or None, (default=None)
            Number of components to keep. If n_components `` is ``None, then n_components is set to min(n_samples, n_features).
        whiten bool, optional
            When True (False by default) the components_ vectors are divided by n_samples times components_ to ensure uncorrelated outputs with unit component-wise variances.
        copy bool, (default=True)
            If False, X will be overwritten. copy=False can be used to save memory but is unsafe for general use.
        batch_sizeint or None, (default=None)
            The number of samples to use for each batch. Only used when calling fit. If batch_size is None, then batch_size is inferred from the data and set to 5 * n_features, to provide a balance between approximation accuracy and memory consumption.
    Attributes
        components_array, shape (n_components, n_features)
            Components with maximum variance.
        explained_variance_array, shape (n_components,)
            Variance explained by each of the selected components.
        explained_variance_ratio_array, shape (n_components,)
            Percentage of variance explained by each of the selected components. If all components are stored, the sum of explained variances is equal to 1.0.
        singular_values_array, shape (n_components,)
            The singular values corresponding to each of the selected components. The singular values are equal to the 2-norms of the n_components variables in the lower-dimensional space.
        mean_array, shape (n_features,)
            Per-feature empirical mean, aggregate over calls to partial_fit.
        var_array, shape (n_features,)
            Per-feature empirical variance, aggregate over calls to partial_fit.
        noise_variance_float
            The estimated noise covariance following the Probabilistic PCA model from Tipping and Bishop 1999. See “Pattern Recognition and Machine Learning” by C. Bishop, 12.2.1 p. 574 or http://www.miketipping.com/papers/met-mppca.pdf.
        n_components_int
            The estimated number of components. Relevant when n_components=None.
        n_samples_seen_int
            The number of samples processed by the estimator. Will be reset on new calls to fit, but increments across partial_fit calls.
        batch_size_int
            Inferred batch size from batch_size.
    Examples
        >>> from sklearn.datasets import load_digits
        >>> from sklearn.decomposition import IncrementalPCA
        >>> from scipy import sparse
        >>> X, _ = load_digits(return_X_y=True)
        >>> transformer = IncrementalPCA(n_components=7, batch_size=200)
        >>> # either partially fit on smaller batches of data
        >>> transformer.partial_fit(X[:100, :])
        IncrementalPCA(batch_size=200, n_components=7)
        >>> # or let the fit function itself divide the data into batches
        >>> X_sparse = sparse.csr_matrix(X)
        >>> X_transformed = transformer.fit_transform(X_sparse)
        >>> X_transformed.shape
        (1797, 7)
    """

方法参考PCA

2.5.1.3. 基于随机化SVD的PCA 通过丢弃具有较低奇异值的奇异向量的分量,将数据降维到低维空间并保留大部分方差信息是非常有意义的。

例如,如果我们使用64x64像素的灰度级图像进行人脸识别,数据的维数为4096,在这样大的数据上训练含RBF内核的支持向量机是很慢的。此外我们知道这些数据固有维度远低于4096,因为人脸的所有照片都看起来有点相似。样本位于较低维度的流体上(例如约200维)。 PCA算法可以用于线性变换数据,同时降低维数并同时保留大部分可描述的方差信息。

在这种情况下,使用可选参数 svd_solver='randomized' 的 PCA 是非常有用的。既然我们将要丢弃大部分奇异值,那么仅仅就实际转换中所需的奇异向量进行计算就可以使得 PCA 计算过程变得异常有效。

例如:以下显示了来自 Olivetti 数据集的 16 个样本肖像(以 0.0 为中心)。右侧是前 16 个奇异向量重画的肖像。因为我们只需要使用大小为 n_{samples} = 400 和 n_{features} = 64 \times 64 = 4096 的数据集的前 16 个奇异向量, 使得计算时间小于 1 秒。

orig_img pca_img

注意:使用可选参数 svd_solver='randomized' ,在 PCA 中我们还需要给出输入低维空间大小 n_components 。

我们注意到, 如果 n_{\max} = \max(n_{\mathrm{samples}}, n_{\mathrm{features}}) 且 n_{\min} = \min(n_{\mathrm{samples}}, n_{\mathrm{features}}), 对于PCA中的实现,随机 PCA 的时间复杂度是:O(n_{\max}^2 \cdot n_{\mathrm{components}}), 而不是 O(n_{\max}^2 \cdot n_{\min}) 。

就内部实现的方法而言, 随机 PCA 的内存占用量和 2 \cdot n_{\max} \cdot n_{\mathrm{components}}, 而不是 n_{\max}\cdot n_{\min} 成正比。

注意:选择参数 svd_solver='randomized' 的 PCA 的 inverse_transform 的实现, 并不是对应 transform 的逆变换(即使 参数设置为默认的 whiten=False)

示例:

Faces recognition example using eigenfaces and SVMs Faces dataset decompositions 参考资料:

“Finding structure with randomness: Stochastic algorithms for constructing approximate matrix decompositions” Halko, et al., 2009 2.5.1.4. 核 PCA KernelPCA 是 PCA 的扩展,通过使用核方法实现非线性降维(dimensionality reduction) (参阅 成对的矩阵, 类别和核函数)。 它具有许多应用,包括去噪, 压缩和结构化预测( structured prediction ) (kernel dependency estimation(内核依赖估计))。 KernelPCA 支持 transform 和 inverse_transform 。

http://sklearn.apachecn.org/cn/0.19.0/_images/sphx_glr_plot_kernel_pca_0011.png

示例:

Kernel PCA 2.5.1.5. 稀疏主成分分析 ( SparsePCA 和 MiniBatchSparsePCA ) SparsePCA 是 PCA 的一个变体,目的是提取能最大程度得重建数据的稀疏分量集合。

小批量稀疏 PCA ( MiniBatchSparsePCA ) 是一个 SparsePCA 的变体,它速度更快但准确度有所降低。对于给定的迭代次数,通过迭代该组特征的小块来达到速度的增加。

Principal component analysis(主成分分析) (PCA) 的缺点在于,通过该方法提取的成分具有独占的密度表达式,即当表示为原始变量的线性组合时,它们具有非零系数,使之难以解释。在许多情况下,真正的基础分量可以被更自然地想象为稀疏向量; 例如在面部识别中,每个分量可能自然地映射到面部的某个部分。

稀疏的主成分产生更节约、可解释的表达式,明确强调了样本之间的差异性来自哪些原始特征。

以下示例说明了使用稀疏 PCA 提取 Olivetti 人脸数据集中的 16 个分量。可以看出正则化项产生了许多零。此外,数据的自然结构导致了非零系数垂直相邻 (vertically adjacent)。该模型并不具备纯数学意义的执行: 每个分量都是一个向量 h \in \mathbf{R}^{4096}, 除非人性化地的可视化为 64x64 像素的图像,否则没有垂直相邻性的概念。下面显示的分量看起来局部化(appear local)是数据的内在结构的影响,这种局部模式使重建误差最小化。有一种考虑到邻接性和不同结构类型的导致稀疏的规范(sparsity-inducing norms),参见 [Jen09] 对这种方法进行了解。有关如何使用稀疏 PCA 的更多详细信息,请参阅下面的示例部分。

pca_img spca_img

请注意,有多种不同的计算稀疏PCA 问题的公式。 这里使用的方法基于 [Mrl09] 。对应优化问题的解决是一个带有惩罚项(L1范数的) \ell_1 的 PCA 问题(dictionary learning(字典学习)):

(U^, V^) = \underset{U, V}{\operatorname{arg\,min\,}} & \frac{1}{2}||X-UV||2^2+\alpha||V||_1 \\text{subject to\,} & ||U_k||_2 = 1 \text{ for all }0 \leq k < n{components}

稀疏推导(sparsity-inducing) \ell_1 规范也可以当训练样本很少时,避免从噪声中拟合分量。可以通过超参数 alpha 来调整惩罚程度(或称稀疏度)。值较小会导致温和的正则化因式分解,而较大的值将许多系数缩小到零。

注意 虽然本着在线算法的精神, MiniBatchSparsePCA 类不实现 partial_fit , 因为在线算法是以特征为导向,而不是以样本为导向。

示例:

Faces dataset decompositions 参考资料:

[Mrl09] “Online Dictionary Learning for Sparse Coding” J. Mairal, F. Bach, J. Ponce, G. Sapiro, 2009 [Jen09] “Structured Sparse Principal Component Analysis” R. Jenatton, G. Obozinski, F. Bach, 2009