网站开发需要资质吗,青海省住房和城乡建设局网站首页,抚顺市建设局网站,wordpress首行缩进PCA#xff08;Principal Component Analysis#xff09; 是一种常见的数据分析方式#xff0c;常用于高维数据的降维#xff0c;可用于提取数据的主要特征分量。
最大可分性
基向量乘原始矩阵会将矩阵映射到这个基向量空间中#xff0c;如果基的数量少于向量本身的维数…PCAPrincipal Component Analysis 是一种常见的数据分析方式常用于高维数据的降维可用于提取数据的主要特征分量。
最大可分性
基向量乘原始矩阵会将矩阵映射到这个基向量空间中如果基的数量少于向量本身的维数则可以达到降维的效果。
如何选取基
希望投影后的投影值尽可能分散因为如果重叠就会有样本消失。当然这个也可以从熵的角度进行理解熵越大所含信息越多。
数据的分散程度可以用方差来表示所以需要将方差最大化
一维是方差而对于高维数据我们用协方差进行约束协方差可以表示两个变量的相关性。为了让两个变量尽可能表示更多的原始信息我们希望它们之间不存在线性相关性。
协方差公式
为了方便处理我们将每个变量的均值都化为 0 当样本数较大时不必在意其是 m 还是 m-1为了方便计算我们分母取 m。
当协方差为0时表示两个变量线性不相关为了让协方差为0我们选择第二个基时只能在与第一个基正交的方向上进行选择因此最终选择的两个方向一定是正交的。
至此我们得到了降维问题的优化目标将一组 N 维向量降为 K 维其目标是选择 K 个单位正交基使得原始数据变换到这组基上后各变量两两间协方差为 0而变量方差则尽可能大在正交的约束下取最大的 K 个方差 根据我们的优化条件我们需要将除对角线外的其它元素化为 0并且在对角线上将元素按大小从上到下排列变量方差尽可能大
设原始数据矩阵 X 对应的协方差矩阵为 C而 P 是一组基按行组成的矩阵设 YPX则 Y 为 X 对 P 做基变换后的数据。设 Y 的协方差矩阵为 D我们推导一下 D 与 C 的关系 所以P为所要求的基矩阵
求解步骤
总结一下 PCA 的算法步骤
设有 m 条 n 维数据。
将原始数据按列组成 n 行 m 列矩阵 X将 X 的每一行进行零均值化即减去这一行的均值求出协方差矩阵 求出协方差矩阵的特征值及对应的特征向量将特征向量按对应特征值大小从上到下按行排列成矩阵取前 k 行组成矩阵 PYPX即为降维到 k 维后的数据
PCA与SVD本质一样
SVD
SVD的基矩阵是和的特征值分解的特征向量按列组成的正交矩阵左奇异矩阵V右奇异矩阵UPCA只与SVD的右奇异向量的压缩效果相同
为什么分左奇异右奇异矩阵
因为SVD所求的矩阵不是方阵协方差矩阵不一样
当矩阵为方阵是SVD等价于PCA
当则
SVD与PCA等价所以PCA问题可以转化为SVD问题求解那转化为SVD问题有什么好处
有三点
一般 X 的维度很高 的计算量很大方阵的特征值分解计算效率不高SVD除了特征值分解这种求解方式外还有更高效且更准确的迭代求解法避免了 的计算
其实PCA只与SVD的右奇异向量的压缩效果相同。
如果取V的前 k 行作为变换矩阵 则 起到压缩行即降维的效果如果取 U的前 k行作为变换矩阵 则 起到压缩列即去除冗余样本的效果。 from __future__ import print_function
from numpy import *
import matplotlib.pyplot as plt
print(__doc__)def loadDataSet(fileName,delim\t):fropen(fileName)stringArr[line.strip().split(delim) for line in fr.readlines()]datArr[map(float,line) for line in stringArr]return mat(datArr)def pca(dataMat,topNfeat9999999):pcaArgs:dataMat 原数据集矩阵topNfeat 应用的N个特征Returns:lowDDataMat 降维后数据集reconMat 新的数据集空间# 计算每一列的均值meanValsmean(dataMat,axis0)# print meanVals, meanVals# 每个向量同时都减去 均值meanRemoveddataMat-meanVals# print meanRemoved, meanRemoved# cov协方差[(x1-x均值)*(y1-y均值)(x2-x均值)*(y2-y均值)...(xn-x均值)*(yn-y均值)]/(n-1)方差: 一维度量两个随机变量关系的统计量协方差: 二维度量各个维度偏离其均值的程度协方差矩阵: 多维度量各个维度偏离其均值的程度当 cov(X, Y)0时表明X与Y正相关(X越大Y也越大X越小Y也越小。这种情况我们称为“正相关”。)当 cov(X, Y)0时表明X与Y负相关当 cov(X, Y)0时表明X与Y不相关。covMatcov(meanRemoved,rowvar0)# eigVals为特征值 eigVects为特征向量eigVals,eigVectslinalg.eig(mat(covMat))# print eigVals, eigVals# print eigVects, eigVects# 对特征值进行从小到大的排序返回从小到大的index序号# 特征值的逆序就可以得到topNfeat个最大的特征向量 x np.array([3, 1, 2]) np.argsort(x)array([1, 2, 0]) # index,1 1; index,2 2; index,0 3 y np.argsort(x) y[::-1]array([0, 2, 1]) y[:-3:-1]array([0, 2]) # 取出 -1, -2 y[:-6:-1]array([0, 2, 1])eigValIndargsort(eigVals)# print eigValInd1, eigValInd# -1表示倒序返回topN的特征值[-1 到 -(topNfeat1) 但是不包括-(topNfeat1)本身的倒叙]eigValIndeigValInd[:-(topNfeat1):-1]# print eigValInd2, eigValInd# 重组 eigVects 最大到最小redEigVectseigVects[:,eigValInd]# print redEigVects, redEigVects.T# 将数据转换到新空间# print ---, shape(meanRemoved), shape(redEigVects)lowDDataMatmeanRemoved*redEigVectsreconMat(lowDDataMat*redEigVects.T)meanVals# print lowDDataMat, lowDDataMat# print reconMat, reconMatreturn lowDDataMat,reconMatdef replaceNanWithMean():datMatloadDataSet(data/13.PCA/secom.data, )numFeatshape(datMat)[1]for i in range(numFeat):# 对value不为NaN的求均值# .A 返回矩阵基于的数组meanValmean(datMat[nonzero(~isnan(datMat[:,i].A))[0],i])# 将value为NaN的值赋值为均值datMat[nonzero(isnan(datMat[:,i].A))[0],i]meanValreturn datMatdef show_picture(dataMat,reconMat):figplt.figure()axfig.add_subplot(111)ax.scatter(dataMat[:,0].flatten().A[0],dataMat[:,1].flatten().A[0],marker^,s90,cgreen)ax.scatter(reconMat[:, 0].flatten().A[0], reconMat[:, 1].flatten().A[0], markerv, s50, cred)plt.show()def analyse_data(dataMat):meanValsmean(dataMat,axis0)meanRemoveddataMat-meanValscovMatcov(meanRemoved,rowvar0)eigvals,eigVectslinalg.eig(mat(covMat))eigValIndargsort(eigvals)topNfeat20eigValIndeigValInd[:-(topNfeat1):-1]cov_all_scorefloat(sum(eigvals))sum_cov_score0for i in range(0,len(eigValInd)):line_cov_scorefloat(eigvals[eigValInd[i]])sum_cov_scoreline_cov_score我们发现其中有超过20%的特征值都是0。这就意味着这些特征都是其他特征的副本也就是说它们可以通过其他特征来表示而本身并没有提供额外的信息。最前面15个值的数量级大于10^5实际上那以后的值都变得非常小。这就相当于告诉我们只有部分重要特征重要特征的数目也很快就会下降。最后我们可能会注意到有一些小的负值他们主要源自数值误差应该四舍五入成0.print(主成分: %s, 方差占比: %s%%, 累积方差占比: %s%% % (format(i1, 2.0f), format(line_cov_score/cov_all_score*100, 4.2f), format(sum_cov_score/cov_all_score*100, 4.1f)))
import numpy as np
xnp.array([3,1,2])
print(np.argsort(x))
y np.argsort(x)
print(y)
print(y[::-1])
print(y[1:3:1])
# 加载数据并转化数据类型为float
dataMatloadDataSet(data/13.PCA/testSet.txt)
# 只需要1个特征向量
lowDmat,reconMatpca(dataMat,1)
# 只需要2个特征向量和原始数据一致没任何变化
lowDmat,reconMatpca(dataMat,2)
print(shape(lowDmat))
show_picture(dataMat,reconMat)
show_picture(dataMat,lowDmat)
show_picture(reconMat,lowDmat)
# 利用PCA对半导体制造数据降维
dataMat replaceNanWithMean()
print(shape(dataMat))
# 分析数据
analyse_data(dataMat)
lowDmat, reconMat pca(dataMat, 20)
print (shape(lowDmat))
show_picture(dataMat,reconMat)
show_picture(dataMat,lowDmat)
show_picture(reconMat,lowDmat) SVD
应用场景
1.信息检索-隐性语义检索Latent Semantic Indexing, LSI或 隐形语义分析Latent Semantic Analysis, LSA
2.推荐系统
利用 SVD 从数据中构建一个主题空间。再在该空间下计算其相似度。(从高维-低维空间的转化在低维空间来计算相似度SVD 提升了推荐系统的效率。
3.图像压缩
将图像矩阵进行奇异值分解然后存储
例如: 32*321024 32*22*132*2130(2*1表示去掉了除对角线的0), 几乎获得了10倍的压缩比。
推荐系统
基于物品的相似度和基于用户的相似度: 物品比较少则选择物品相似度用户比较少则选择用户相似度。
基于物品的相似度: 计算物品之间的距离。【耗时会随物品数量的增加而增加】由于物品A和物品C 相似度(相关度)很高所以给买A的人推荐C。
基于用户的相似度: 计算用户之间的距离。【耗时会随用户数量的增加而增加】由于用户A和用户C 相似度(相关度)很高所以A和C是兴趣相投的人对于C买的物品就会推荐给A。
相似度计算
inA, inB 对应的是 列向量
欧氏距离: 指在m维空间中两个点之间的真实距离或者向量的自然长度即该点到原点的距离。二维或三维中的欧氏距离就是两点之间的实际距离。 相似度 1/(1欧式距离)相似度 1.0/(1.0 la.norm(inA - inB))物品对越相似它们的相似度值就越大。皮尔逊相关系数: 度量的是两个向量之间的相似度。 相似度 0.5 0.5*corrcoef() 【皮尔逊相关系数的取值范围从 -1 到 1通过函数0.5 0.5*corrcoef()这个函数计算把值归一化到0到1之间】相似度 0.5 0.5 * corrcoef(inA, inB, rowvar 0)[0][1]相对欧氏距离的优势: 它对用户评级的量级并不敏感。余弦相似度: 计算的是两个向量夹角的余弦值。 余弦值 (A·B)/(||A||·||B||) 【余弦值的取值范围也在-1到1之间】相似度 0.5 0.5*余弦值相似度 0.5 0.5*( float(inA.T*inB) / la.norm(inA)*la.norm(inB))如果夹角为90度则相似度为0如果两个向量的方向相同则相似度为1.0。
推荐系统的原理
推荐系统的工作过程: 给定一个用户系统会为此用户返回N个最好的推荐菜。实现流程大致如下: 寻找用户没有评级的菜肴即在用户-物品矩阵中的0值。在用户没有评级的所有物品中对每个物品预计一个可能的评级分数。这就是说: 我们认为用户可能会对物品的打分这就是相似度计算的初衷。对这些物品的评分从高到低进行排序返回前N个物品。
#SVD
from __future__ import print_function
from numpy import linalg as la
from numpy import *def loadExData3():# 利用SVD提高推荐效果菜肴矩阵return[[2, 0, 0, 4, 4, 0, 0, 0, 0, 0, 0],[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 5],[0, 0, 0, 0, 0, 0, 0, 1, 0, 4, 0],[3, 3, 4, 0, 3, 0, 0, 2, 2, 0, 0],[5, 5, 5, 0, 0, 0, 0, 0, 0, 0, 0],[0, 0, 0, 0, 0, 0, 5, 0, 0, 5, 0],[4, 0, 4, 0, 0, 0, 0, 0, 0, 0, 5],[0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 4],[0, 0, 0, 0, 0, 0, 5, 0, 0, 5, 0],[0, 0, 0, 3, 0, 0, 0, 0, 4, 5, 0],[1, 1, 2, 1, 1, 2, 1, 0, 4, 5, 0]]def loadExData2():# 书上代码给的示例矩阵return[[0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 5],[0, 0, 0, 3, 0, 4, 0, 0, 0, 0, 3],[0, 0, 0, 0, 4, 0, 0, 1, 0, 4, 0],[3, 3, 4, 0, 0, 0, 0, 2, 2, 0, 0],[5, 4, 5, 0, 0, 0, 0, 5, 5, 0, 0],[0, 0, 0, 0, 5, 0, 1, 0, 0, 5, 0],[4, 3, 4, 0, 0, 0, 0, 5, 5, 0, 1],[0, 0, 0, 4, 0, 4, 0, 0, 0, 0, 4],[0, 0, 0, 2, 0, 2, 5, 0, 0, 1, 2],[0, 0, 0, 0, 5, 0, 0, 0, 0, 4, 0],[1, 0, 0, 0, 0, 0, 0, 1, 2, 0, 0]]
def loadExData():# 推荐引擎示例矩阵return[[4, 4, 0, 2, 2],[4, 0, 0, 3, 3],[4, 0, 0, 1, 1],[1, 1, 1, 2, 0],[2, 2, 2, 0, 0],[1, 1, 1, 0, 0],[5, 5, 5, 0, 0]]# # 原矩阵# return[[1, 1, 1, 0, 0],# [2, 2, 2, 0, 0],# [1, 1, 1, 0, 0],# [5, 5, 5, 0, 0],# [1, 1, 0, 2, 2],# [0, 0, 0, 3, 3],# [0, 0, 0, 1, 1]]# 原矩阵return[[0, -1.6, 0.6],[0, 1.2, 0.8],[0, 0, 0],[0, 0, 0]]# 相似度计算假定inA和inB 都是列向量
# 基于欧氏距离
def ecludSim(inA,inB):return 1.0/(1.0la.norm(inA-inB))# pearsSim()函数会检查是否存在3个或更多的点。
# corrcoef直接计算皮尔逊相关系数范围[-1, 1]归一化后[0, 1]
def pearsSim(inA,inB):# 如果不存在该函数返回1.0此时两个向量完全相关。if len(inA)3:return 1.0return 0.50.5*corrcoef(inA,inB,rowvar0)[0][1]# 计算余弦相似度如果夹角为90度相似度为0如果两个向量的方向相同相似度为1.0
def cosSim(inA,inB):numfloat(inA.T*inB)denomla.norm(inA)*la.norm(inB)return 0.50.5*(num/denom)# 基于物品相似度的推荐引擎
def standEst(dataMat,user,simMeas,item):standEst(计算某用户未评分物品中以对该物品和其他物品评分的用户的物品相似度然后进行综合评分)Args:dataMat 训练数据集user 用户编号simMeas 相似度计算方法item 未评分的物品编号Returns:ratSimTotal/simTotal 评分05之间的值# 得到数据集中的物品数目nshape(dataMat)[1]# 初始化两个评分值simTotal0.0ratSimTotal0.0# 遍历行中的每个物品对用户评过分的物品进行遍历并将它与其他物品进行比较for j in range(n):userRatingdataMat[user,j]# 如果某个物品的评分值为0则跳过这个物品if userRating0:continue# 寻找两个用户都评级的物品# 变量 overLap 给出的是两个物品当中已经被评分的那个元素的索引ID# logical_and 计算x1和x2元素的真值。overLapnonzero(logical_and(dataMat[:,item].A0,dataMat[:,j].A0))[0]# 如果相似度为0则两着没有任何重合元素终止本次循环if len(overLap)0:similarity0# 如果存在重合的物品则基于这些重合物重新计算相似度。else:similaritysimMeas(dataMat[overLap,item],dataMat[overLap,j])# print the %d and %d similarity is : %f(iten,j,similarity)# 相似度会不断累加每次计算时还考虑相似度和当前用户评分的乘积# similarity 用户相似度 userRating 用户评分simTotalsimilarityratSimTotalsimilarity*userRatingif simTotal0:return 0# 通过除以所有的评分总和对上述相似度评分的乘积进行归一化使得最后评分在0~5之间这些评分用来对预测值进行排序else:return ratSimTotal/simTotal# 基于SVD的评分估计
# 在recommend() 中这个函数用于替换对standEst()的调用该函数对给定用户给定物品构建了一个评分估计值
def svdEst(dataMat,user,simMeas,item):svdEst( )Args:dataMat 训练数据集user 用户编号simMeas 相似度计算方法item 未评分的物品编号Returns:ratSimTotal/simTotal 评分05之间的值# 物品数目nshape(dataMat)[1]# 对数据集进行SVD分解simTotal0.0ratSimTotal0.0# 奇异值分解# 在SVD分解之后我们只利用包含了90%能量值的奇异值这些奇异值会以NumPy数组的形式得以保存U,Sigma,VTla.svd(dataMat)# # 分析 Sigma 的长度取值# analyse_data(Sigma, 20)# 如果要进行矩阵运算就必须要用这些奇异值构建出一个对角矩阵Sig4mat(eye(4)*Sigma[:4])# 利用U矩阵将物品转换到低维空间中构建转换后的物品(物品4个主要的特征)xformedItems dataMat.T * U[:, :4] * Sig4.IxformedItems1 U[:, :4].T *dataMatprint(dataMat, shape(dataMat))print(U[:, :4], shape(U[:, :4]))print(Sig4.I, shape(Sig4.I))print(VT[:4, :], (VT[:4, :]))print(xformedItems, (xformedItems))print(xformedItems1, (xformedItems1))# 对于给定的用户for循环在用户对应行的元素上进行遍历# 这和standEst()函数中的for循环的目的一样只不过这里的相似度计算时在低维空间下进行的。for j in range(n):userRating dataMat[user, j]if userRating 0 or j item:continue# 相似度的计算方法也会作为一个参数传递给该函数similarity simMeas(xformedItems[item, :].T, xformedItems[j, :].T)# for 循环中加入了一条print语句以便了解相似度计算的进展情况。如果觉得累赘可以去掉print(the %d and %d similarity is: %f % (item, j, similarity))# 对相似度不断累加求和simTotal similarity# 对相似度及对应评分值的乘积求和ratSimTotal similarity * userRatingif simTotal 0:return 0else:# 计算估计评分return ratSimTotal/simTotal# recommend()函数就是推荐引擎它默认调用standEst()函数产生了最高的N个推荐结果。
# 如果不指定N的大小则默认值为3。该函数另外的参数还包括相似度计算方法和估计方法
def recommend(dataMat,user,N3,simMeascosSim,estMethodstandEst):svdEst( )Args:dataMat 训练数据集user 用户编号simMeas 相似度计算方法estMethod 使用的推荐算法Returns:返回最终 N 个推荐结果# 寻找未评级的物品# 对给定的用户建立一个未评分的物品列表unratedItemsnonzero(dataMat[user,:].A0)[1]# 如果不存在未评分物品那么就退出函数if len(unratedItems)0:return you rated everything# 物品的编号和评分值itemScores[]# 在未评分物品上进行循环for item in unratedItems:# 获取 item 该物品的评分estimatedScoreestMethod(dataMat,user,simMeas,item)itemScores.append((item,estimatedScore))# 按照评分得分 进行逆排序获取前N个未评级物品进行推荐return sorted(itemScores,keylambda jj:jj[1],reverseTrue)[:N]def analyse_data(Sigma, loopNum20):analyse_data(分析 Sigma 的长度取值)Args:Sigma Sigma的值loopNum 循环次数# 总方差的集合总能量值Sig2 Sigma**2SigmaSum sum(Sig2)for i in range(loopNum):SigmaI sum(Sig2[:i1])根据自己的业务情况就行处理设置对应的 Singma 次数通常保留矩阵 80% 90% 的能量就可以得到重要的特征并取出噪声。print(主成分: %s, 方差占比: %s%% % (format(i1, 2.0f), format(SigmaI/SigmaSum*100, 4.2f)))# 图像压缩函数
# 加载并转换数据
def imgLoadData(filename):myl[]# 打开文本文件并从文件以数组方式读入字符for line in open(filename).readlines():newRow[]for i in range(32):newRow.append(int(line[i]))myl.append(newRow)# 矩阵调入后就可以在屏幕上输出该矩阵myMat mat(myl)return myMat# 打印矩阵
def printMat(inMat, thresh0.8):# 由于矩阵保护了浮点数因此定义浅色和深色遍历所有矩阵元素当元素大于阀值时打印1否则打印0for i in range(32):for k in range(32):if float(inMat[i, k]) thresh:print(1, end )else:print(0, end )print()# 实现图像压缩允许基于任意给定的奇异值数目来重构图像
def imgCompress(numSV3, thresh0.8):imgCompress( )Args:numSV Sigma长度 thresh 判断的阈值# 构建一个列表myMatimgLoadData(data/14.SVD/0_5.txt)print(****original matrix****)# 对原始图像进行SVD分解并重构图像eprintMat(myMat, thresh)# 通过Sigma 重新构成SigRecom来实现# Sigma是一个对角矩阵因此需要建立一个全0矩阵然后将前面的那些奇异值填充到对角线上。U, Sigma, VT la.svd(myMat)# SigRecon mat(zeros((numSV, numSV)))# for k in range(numSV):# SigRecon[k, k] Sigma[k]# 分析插入的 Sigma 长度analyse_data(Sigma, 20)SigRecon mat(eye(numSV) * Sigma[: numSV])reconMat U[:, :numSV] * SigRecon * VT[:numSV, :]print(****reconstructed matrix using %d singular values ***** % numSV)printMat(reconMat, thresh)
# 对矩阵进行SVD分解(用python实现SVD)
Data loadExData()
print (Data:, Data)
U, Sigma, VT linalg.svd(Data)
# 打印Sigma的结果因为前3个数值比其他的值大了很多为9.72140007e005.29397912e006.84226362e-01
# 后两个值比较小每台机器输出结果可能有不同可以将这两个值去掉
print (U:, U)
print (Sigma, Sigma)
print (VT:, VT)
print (VT:, VT.T)
# 重构一个3x3的矩阵Sig3
Sig3 mat([[Sigma[0], 0, 0], [0, Sigma[1], 0], [0, 0, Sigma[2]]])
print (U[:, :3] * Sig3 * VT[:3, :])
# 计算欧氏距离
myMat mat(loadExData())
# print myMat
print (ecludSim(myMat[:, 0], myMat[:, 1]))
print (ecludSim(myMat[:, 0], myMat[:, 0]))# 计算余弦相似度
print (cosSim(myMat[:, 0], myMat[:, 1]))
print (cosSim(myMat[:, 0], myMat[:, 0]))# 计算皮尔逊相关系数
print (pearsSim(myMat[:, 0], myMat[:, 1]))
print (pearsSim(myMat[:, 0], myMat[:, 0]))# 计算相似度的方法
myMat mat(loadExData3())
# print myMat
# 计算相似度的第一种方式
print(recommend(myMat, 1, estMethodsvdEst))
# 计算相似度的第二种方式
print(recommend(myMat, 1, estMethodsvdEst, simMeaspearsSim))
# 默认推荐菜馆菜肴推荐示例
print(recommend(myMat, 2))# 利用SVD提高推荐效果
U, Sigma, VT la.svd(mat(loadExData2()))
print (Sigma) # 计算矩阵的SVD来了解其需要多少维的特征
Sig2 Sigma**2 # 计算需要多少个奇异值能达到总能量的90%
print (sum(Sig2)) # 计算总能量
print (sum(Sig2) * 0.9) # 计算总能量的90%
print (sum(Sig2[: 2])) # 计算前两个元素所包含的能量
print (sum(Sig2[: 3])) # 两个元素的能量值小于总能量的90%于是计算前三个元素所包含的能量
# 该值高于总能量的90%这就可以了# 压缩图片
imgCompress(3)
