pca 源码

#!/usr/python/bin
#-*- coding:utf-8 -*-

from numpy import *
import matplotlib
import matplotlib.pyplot as plt

def loadDataSet(fileName, delim='\t'):
	fr = open(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, topNfeat = 9999999):			#输入原始数据矩阵，返回降维之后的数据矩阵
	meanVals = mean(dataMat, axis = 0)
	meanRemoved = dataMat-meanVals				#去除平均值，因为协方差矩阵的定义就要每个元素减去其平均值
	covMat = cov(meanRemoved, rowvar = 0)		#计算协方差矩阵，协方差矩阵每一维的大小就是原数据中特征的数量
	eigVals, eigVects = linalg.eig(mat(covMat))	#计算协方差矩阵的特征值和特征向量
	eigValInd = argsort(eigVals)				#返回的特征值本来就是有序的(反正是取最大的几个特征值对应的特征向量)
	eigValInd = eigValInd[:-(topNfeat+1):-1]	#返回倒数第topNfeat个特征值对应的下标(反正是取最大的几个特征值对应的特征向量)
	redEigVects = eigVects[:,eigValInd]			#根据特征值的排序结果，选出对应重要程度排名前几位的特征向量
	lowDDataMat = meanRemoved * redEigVects		#计算降维之后的数据
	reconMat = (lowDDataMat * redEigVects.T)+meanVals 	#将降维之后的数据映射到原空间
	return lowDDataMat, reconMat	#降维之后的数据和映射到原空间的数据

def pic(dataMat,reconMat):
	fig = plt.figure()
	ax = fig.add_subplot(111)
	ax.scatter(dataMat[:,0].flatten().A[0], dataMat[:,1].flatten().A[0],marker = '^',s = 90)
	ax.scatter(reconMat[:,0].flatten().A[0], reconMat[:,1].flatten().A[0],marker = 'o',s = 50,c = 'red')
	fig.show()

def replaceNanWithMean():		#使用每个特征的平均值代替特征中的缺失值，返回填补完整的数据矩阵
	datMat = loadDataSet('secom.data',' ')
	numFeat = shape(datMat)[1]
	for i in range(numFeat):		#因为要计算每个特征的平均值，所以要取总的特征数量
		meanVal = mean(datMat[nonzero(~isnan(datMat[:,i].A))[0],i])		#计算每个特征的平均值（出去NAN的值）
		datMat[nonzero(isnan(datMat[:,i].A))[0],i] = meanVal 		#将NAN的值替换为相应属性的
	return datMat

协方差矩阵的示意图

PCA的主要步骤：

1. 去除每个特征的平均值（计算协方差矩阵的需要）

2. 计算协方差矩阵（协方差矩阵的定义）

3. 计算协方差矩阵的特征值和特征向量

4. 根据特征值的大小排序

5. 选取最大的几个特征值对应的特征向量构成特征矩阵

6. 原数据数组与特征矩阵相乘，得到降维之后的数据

7. 可以将降维之后的数据复原到原空间中（乘特征矩阵的逆）与原数组画图做对比（低维度才能看清楚）