使用SVD奇异值分解求解PCA+Python实现
这几天在看有关PCA的博客时,感觉文章中针对如何用SVD解PCA的过程的讲解不是很清晰,自己捋了捋思路,把自己对于SVD解PCA的步骤分享出来
关于PCA的思想和原理,这里不阐述,推荐这篇文章,当初就是看这篇文章入坑PCA的
http://blog.codinglabs.org/articles/pca-tutorial.html
1.直接解PCA
设有n组m维数据,组成 m*n 的矩阵 Xm∗n\Chi_{m*n}Xm∗n , 则 在去中心化后,X\ChiX 的协方差矩阵
C\mathcal{C}C = 1n\frac{1}{n}n1 X\ChiXXT\Chi^{T}XT
目标:将 C\mathcal{C}C 对角化
如何对角化 C\mathcal{C}C:
设P\mathcal{P}P是正交基按行组成的矩阵(最终希望得到的转换矩阵,P\mathcal{P}P的前k\mathcal{k}k行就是要寻找的基,可对X\ChiX降维),则Y\mathcal{Y}Y=P\mathcal{P}PX\ChiX中,Y\mathcal{Y}Y是X\ChiX对P\mathcal{P}P做基变换得到的结果
则Y\mathcal{Y}Y 的协方差矩阵D\mathcal{D}D = 1n\frac{1}{n}n1Y\mathcal{Y}YYT\mathcal{Y}^{T}YT = P\mathcal{P}PC\mathcal{C}CPT\mathcal{P}^{T}PT 应该是一个对角阵。
已知知道一个m∗\mathcal{m}*m∗m\mathcal{m}m的实对称矩阵一定可以找到m\mathcal{m}m个正交特征向量,特征向量组成的矩阵为E\mathcal{E}E
则有 ET\mathcal{E}^{T}ETC\mathcal{C}CE\mathcal{E}E = Λ\LambdaΛ,对比上式,可以得到 P\mathcal{P}P = ET\mathcal{E}^{T}ET
即我们现在的目标是求出C\mathcal{C}C的特征向量矩阵,并将其转置
但是,高维数据的协方差矩阵非常庞大,比如单位样本是一副100*100的uint8图像,则一副图像就要表示为一个 10000 * 1 的列向量,当样本数量少于10000时,协方差矩阵的大小为10000 *10000(800m左右),计算机求解这个矩阵会非常困难,我尝试过,用Matlab求解要用5分钟左右的时间,如果图像更大,比如时 70000 * 1 的情况,则协方差矩阵大小达到45G,求解成为了几乎不可能的事情
这时候我们就要借助一个工具 SVD特征值分解,来曲线求解协方差矩阵的特征值矩阵
2. 用SVD特征值分解来解PCA
这里不阐述SVD的原理,同样我推荐这篇文章
https://mp.weixin.qq.com/s/Dv51K8JETakIKe5dPBAPVg
X\ChiX经过分解后,可以得到
Xm∗n\Chi_{m*n}Xm∗n = Um∗m\mathcal{U}_{m*m}Um∗mΣm∗n\Sigma_{m*n}Σm∗nVn∗nT\mathcal{V}^{T}_{n*n}Vn∗nT
其中Um∗m\mathcal{U}_{m*m}Um∗m是X\ChiXXT\Chi^{T}XT的特征向量按列组成的矩阵
Vn∗n\mathcal{V}_{n*n}Vn∗n是XT\Chi^{T}XTX\ChiX的特征向量按列组成的矩阵
Σm∗n\Sigma_{m*n}Σm∗n的对角线上是奇异值σ\sigmaσ,奇异值与特征值的关系为σ2\sigma^{2}σ2 = λ\lambdaλ
按以上的性质,可以得到 P\mathcal{P}P = UT\mathcal{U}^{T}UT
所以我们现在的目标是求出U\mathcal{U}U
SVD的分解思想是
Xm∗n\Chi_{m*n}Xm∗n≈\approx≈Um∗r\mathcal{U}_{m*r}Um∗rΣr∗r\Sigma_{r*r}Σr∗rVr∗nT\mathcal{V}^{T}_{r*n}Vr∗nT
我们现在将Xm∗n\Chi_{m*n}Xm∗n近似为
Xm∗n\Chi_{m*n}Xm∗n≈\approx≈Um∗n\mathcal{U}_{m*n}Um∗nΣn∗n\Sigma_{n*n}Σn∗nVn∗nT\mathcal{V}^{T}_{n*n}Vn∗nT (对于n <= m)
Vn∗n\mathcal{V}_{n*n}Vn∗n的特征值矩阵为Σn∗n2\Sigma_{n*n}^{2}Σn∗n2
所以我们可以求得
Um∗n\mathcal{U}_{m*n}Um∗n = Xm∗n\Chi_{m*n}Xm∗n(Σn∗n(\Sigma_{n*n}(Σn∗nVn∗n)−1\mathcal{V}_{n*n})^{-1}Vn∗n)−1
从而就可以依据特征值大小将Um∗n\mathcal{U}_{m*n}Um∗n的特征值排序(U\mathcal{U}U和V\mathcal{V}V拥有相同的特征值),就可以降维至Um∗k\mathcal{U}_{m*k}Um∗k,再将其转置就可以得到变换矩阵P\mathcal{P}P
这时可能会有疑问
Um∗m\mathcal{U}_{m*m}Um∗m和Vn∗n\mathcal{V}_{n*n}Vn∗n维数不一样,为什么特征值个数是min{m, n},其实,若直接对前者进行特征值求解,排序后,会发现,前者只有前n个特征值是实数,后面都是虚数
可以看到,在用SVD求解PCA时,我们只需要计算Vn∗n\mathcal{V}_{n*n}Vn∗n的特征值和特征向量,通常情况下,数据样本 n 是远小于数据维数 m 的,因此计算 n*n 矩阵的特征值和特征向量,将可以大大减小计算量,SVD“曲线救国”的目的也就达到了
3.Python实现PCA
在实际进行PCA的求解中,需要根据数据的维度与样本数来进行直接求解PCA与用SVD求解PCA的权衡
import numpy as npclass PCA():"""PCA降维"""def __init__(self):self.X = Noneself.k = Nonedef getTransMat(self, X, k, svd):"""训练得到特征空间中的主成分转换矩阵args:X -- 训练数据 shape=(features, samples)k -- 要降维的维度,若小于0,则自动选择, 选择阈值为特征值总和95%svd -- 是否使用svd求解PCA, true or false 如果训练样本数多于特征数,选择直接解更省内存(faster),反之选择svd奇异值分解更好Returns:transMat -- 转换矩阵 shape=(降维后维度, 降维前维度)"""# 将X的每一维进行零均值化self.X = Xself.k = km = np.mean(self.X, 1)m = np.array([m])self.X = self.X - m.Tif svd:# 用SVD特征值分解解出变换矩阵# 求出 X' * X 的特征值矩阵V,并求出奇异值对角阵Σ(奇异值从大到小排序)C = np.dot(self.X.T, self.X)# 求出 X' * X 的特征值矩阵及对应的特征向量和奇异值矩阵[eigValue, V] = np.linalg.eig(C)eigValue = abs(np.array([eigValue]))sinValue = np.sqrt(eigValue)temp = np.zeros([np.size(sinValue, 1), np.size(sinValue, 1)])for i in range(0, np.size(sinValue, 1)):temp[i, i] = sinValue[0, i]sinValue = temp# 将特征向量归一化V = V / np.linalg.norm(V, axis=0)# 求出变换矩阵U = np.dot(self.X, np.linalg.inv(np.dot(sinValue, V.T)))U = U.Telse:# 直接解PCAC = np.dot(self.X, self.X.T)[eigValue, V] = np.linalg.eig(C)eigValue = np.array([eigValue])# 将特征向量归一化V = V / np.linalg.norm(V, axis=0)# 求出变换矩阵U = V.T# 按特征值大小给变换矩阵行向量重新排序index = np.argsort(-eigValue)eigValue.tolist()[0].reverse()temp = U.copy()for i in range(0, np.size(U, 0)):U[i, :] = temp[index[0][i], :]# 自动选择降维维度self.autoChooseDim(eigValue)# 获得降维变换矩阵transMat = U[0:self.k, :]return np.real(transMat)def autoChooseDim(self, sortV):"""自适应选择降维维度args:sortV -- 已经排好序的特征值序列returns:None"""if self.k <= 0:threshold = sum(sum(sortV)) * 0.95total = 0for i in range(0, len(sortV[0])):total += sortV[0][i]if total > threshold:breakself.k = ielse:passif __name__ == "__main__":pass
使用SVD奇异值分解求解PCA+Python实现相关推荐
- SVD奇异值分解(PCA,LSI)
版权声明: 本文由LeftNotEasy发布于http://leftnoteasy.cnblogs.com, 本文可以被全部的转载或者部分使用,但请注明出处,如果有问题,请联系wheeleast@gm ...
- 【数学和算法】最小二乘法,SVD奇异值分解、LU分解的应用场景
多项式拟合曲线时使用最小二乘法,把问题化简为A*x=B的线性方程组: 然后使用LU矩阵分解算法求解线性方程组A*x=B,具体做法是: 先对矩阵A进行行初等变换得到上三角矩阵U: 再求出下三角矩阵L,就 ...
- svd奇异值分解_传统推荐算法(一)SVD推荐(1)解读奇异值分解
文章目录 写在前面 1. 从几何变换到奇异值分解 2. 代数角度理解奇异值与奇异向量 2.1 从正交基映射推导SVD 2.2 特征值分解求解奇异值和奇异向量 2.2.1 求解过程 2.2.2 推论 2 ...
- SVD奇异值分解(标题重复率过高)
本篇为<深度学习>系列博客的第三篇,该系列博客主要记录深度学习相关知识的学习过程和自己的理解,方便以后查阅. 之前在看slam的时候就遇到过SVD但是没有做系统的学习,现在接触到深度学习, ...
- svd奇异值分解_SVD(奇异值分解)到底在干什么
奇异值分解就是在低维空间中寻找最接近原矩阵 的低维矩阵 ,说白了就是数据降维. 奇异值分解是一种十分重要但又难以理解的矩阵处理技术,据人工智能的大牛吴恩达老师所说,在机器学习中是最重要的分解没有之一的 ...
- svd奇异值分解_Lecture 28 | 奇异值分解
01 Singular Value Decomposition 奇异值分解 奇异值分解指任一mxn的矩阵A都可以分解为一个mxm酉矩阵U乘一个mxn对角阵Σ再乘一个nxn酉矩阵V共轭转置的形式. 下面 ...
- SVD的原理及python实现——正本清源
SVD的原理及python实现 希望这篇文章能达到正本清源的目的.
- 三硬币问题建模及Gibbs采样求解(Python实现)
大纲 三硬币问题建模及Gibbs采样求解(Python实现) 三硬币问题 二项分布与Beta分布 基础概念 二项分布 Beta分布 二项分布-Beta分布 三硬币问题建模过程 Gibbs 采样求解 p ...
- python 字典处理_python numpy求解积分python中的字典操作及字典函数
字典 dict_fruit = {'apple':'苹果','banana':'香蕉','cherry':'樱桃','avocado':'牛油果','watermelon':'西瓜'} 字典的操作 W ...
- 大M法的python编程求解和python包求解
大M法的python编程求解和python包求解 一.大M算法的求解步骤讲解 二.python编程求解 三.利用python包scipy的优化包optimize 四.用excel求解 五.分析结果 一 ...
最新文章
- Entity Framework简介
- pfx证书密码怎么查询_2019成考成绩查询时间通知!忘记账号和密码怎么办?
- 2.IDA-数据显示窗口(反汇编窗口、函数窗口、十六进制窗口)
- aws emr 大数据分析_DataOps —使用AWS Lambda和Amazon EMR的全自动,低成本数据管道
- hashtable和hashmap的区别?
- Python 内置模块之 asyncio(异步iO)
- (王道408考研操作系统)第三章内存管理-第二节3:页面置换算法2
- 【Linux入门学习之】grep命令详解
- MySQL绿色版mysql-5.7.17-winx64简洁安装教程
- 成人教育考试报名照片的尺寸是多少?大一寸照片怎么做?
- Maze环境以及DQN的实现
- access里面的表达式运用_Access 如何使用表达式生成器
- 2018-01-05-医药行业的IT革命探讨
- LibVLC —— 常用对象解析
- 运行阶段及面向对象技巧
- 服务器内存不足导致程序(tomcat)崩溃
- 活久见!这么好的图文电子书制作工具我竟然才发现
- 备忘3:爬取博主热门信息以及所有热门微博的点赞的用户信息
- C# 科学计数法的转换
- implode() 函数