主成分分析的Python实现

       之前写的博客中提介绍了PCA 的原理,现在这里用Python实现一下吧~所用IDE为Spyder,Python版本3.5。这篇文章中用Mnist作为实验数据集,为了省事,并没有从官网中下载原始的Mnist数据,而是从kaggle上下载的。。上面的train数据和test数据(严格来说应该是验证集,并没有给出标签,而是让我们自己预测的)分别以csv格式给出,不过原始数据有点多,PC跑一些吃不消,所以只用了train中的一部分数据进行试验~首先是原始数据的读取以及初始化工作,train_num、test_num和total_num分别表示训练样本数、测试样本数以及总体样本数,这里的demo我用了其中10000个作为训练样本,10000个作为测试样本。dim为降维后的维数,也就是之前文章提到的投影向量的个数。

"""
reated on Fri May  5 15:46:54 2017

@author: zhizec
"""
import numpy as np
total_num = 20000
train_num = 10000
test_num = 10000
dim = 50
data = np.genfromtxt('D:\\backup_2016107\\train (1).csv',delimiter=',')
number = data.shape[0] - 1
label = data[1:, 0]
train_label = label[0:train_num]
test_label = label[train_num:total_num]
train_samples = data[1:train_num + 1, 1:]
test_samples = data[(train_num + 1):total_num + 1, 1:]

然后是协方差矩阵以及相应特征向量的计算,比较简单:

mean_sample = np.mean(train_samples, axis = 0)
mean_samples = train_samples - np.tile(mean_sample, (train_num,1))
#cal covariance matrix
cov_matrix = np.dot(mean_samples.transpose(), mean_samples)
#cal projection vector
(eigvalue, eigvector) = np.linalg.eig(cov_matrix)
projMatrix = eigvector[:,:dim]
#cal new samples: Y = XW
train_samples_new = np.dot(train_samples, projMatrix)
test_samples_new = np.dot(test_samples, projMatrix)

使用最近邻分类器进行分类,分别计算每个测试样本离 训练样本的距离,以最短距离对应训练样本的标签作为分类的类别,其中index保存了正确分类的个数。在写这一部分代码的时候出现了一些问题,一个就是距离向量的计算,一开始的实现是批量计算欧氏距离,使用tile函数构造出train_num X test_sample的矩阵,然后使用train_samples减去这个距离矩阵得到一个距离矩阵,对每一行求norm就得到了这一行对应训练样本和测试样本的距离,不料跑的时候直接报memory error了,所以只能使用循环分开计算距离了,最后使用argmin求最近的标签。

#NN classifier
index = 0
for i in range(test_num):
    test_sample = test_samples_new[i,:]
    #test_matrix = np.tile(test_samples_new , (10000, 1))
    disctance_vector = np.zeros((train_num,1))
    for j in range(train_num):
        train_sample = train_samples_new[j,:]
        #distance = np.linalg.norm(test_sample - train_sample)
        distance = sum((test_sample - train_sample)**2)
        disctance_vector[j] = distance
    #disctance_matrix = train_samples_new  - test_matrix 
    #disctance_vector = np.linalg.norm(disctance_matrix, axis = 1)
    #sample_label = test_label[i]
    minIndex = np.argmin(disctance_vector)
    print ("actual label: " + repr(train_label[minIndex]) + "expect label: " + repr(test_label[i]))
    if train_label[minIndex] ==  test_label[i]:
        index = index + 1

为了作为比较,给出直接使用最近邻分类器分类的结果,这个时候直接使用原始的训练、测试集合就行了:

index2 = 0 
for i in range(test_num):
    test_sample = test_samples[i,:]
    #test_matrix = np.tile(test_samples_new , (10000, 1))
    disctance_vector = np.zeros((train_num,1))
    for j in range(train_num):
        train_sample = train_samples[j,:]
        #distance = np.linalg.norm(test_sample - train_sample)
        distance = sum((test_sample - train_sample)**2)
        disctance_vector[j] = distance
    #disctance_matrix = train_samples_new  - test_matrix 
    #disctance_vector = np.linalg.norm(disctance_matrix, axis = 1)
    #sample_label = test_label[i]
    minIndex = np.argmin(disctance_vector)
    print ("actual label: " + repr(train_label[minIndex]) + "expect label: " + repr(test_label[i]))
    if train_label[minIndex] ==  test_label[i]:
        index2 = index2 + 1

最后打印出分类识别率:

print (index*1.0/test_num)         
print (index2*1.0/test_num)

运行结果:....actual label: 5.0expect label: 5.0actual label: 6.0expect label: 6.0actual label: 2.0expect label: 2.0actual label: 7.0expect label: 7.0actual label: 3.0expect label: 3.0actual label: 4.0expect label: 4.0actual label: 6.0expect label: 6.0actual label: 1.0expect label: 1.0actual label: 3.0expect label: 5.0actual label: 0.0expect label: 0.0actual label: 4.0expect label: 4.0actual label: 7.0expect label: 7.0actual label: 7.0expect label: 9.00.96480.9524PCA的识别率和降维维数dim有一定的关系,在实际使用中需要通过调节这个值来达到最优值。最后我用kaggle上提供的42000个训练样本对28000个样本进行预测,降到50维,提交之后识别率为0.97514,在训练样本这么大的情况下效果不是特别的高~

# -*- coding: utf-8 -*-
"""
Created on Thu Jun  1 10:14:13 2017

@author: zhizec
"""

total_num = 70000
train_num = 42000
test_num = 28000
dim = 50

import numpy as np
train_data = np.genfromtxt('D:\\Python\\python\\工作相关\\train (1).csv',delimiter=',')
test_data = np.genfromtxt('D:\\Python\\python\\工作相关\\test (1).csv',delimiter=',')
train_label = train_data[1:, 0]

train_samples = train_data[1:,1:]
test_samples = test_data[1:,:]

mean_sample = np.mean(train_samples, axis = 0)
mean_samples = train_samples - np.tile(mean_sample, (train_num,1))
#cal covariance matrix
cov_matrix = np.dot(mean_samples.transpose(), mean_samples)
#cal projection vector
(eigvalue, eigvector) = np.linalg.eig(cov_matrix)
projMatrix = eigvector[:,:dim]
#cal new samples: Y = XW
train_samples_new = np.dot(train_samples, projMatrix)
test_samples_new = np.dot(test_samples, projMatrix)


#NN classifier
test_label = np.zeros((test_num,1))
for i in range(test_num):
    test_sample = test_samples_new[i,:]
    #test_matrix = np.tile(test_samples_new , (10000, 1))
    disctance_vector = np.zeros((train_num,1))
    for j in range(train_num):
        train_sample = train_samples_new[j,:]
        #distance = np.linalg.norm(test_sample - train_sample)
        distance = sum((test_sample - train_sample)**2)
        disctance_vector[j] = distance
    #disctance_matrix = train_samples_new  - test_matrix 
    #disctance_vector = np.linalg.norm(disctance_matrix, axis = 1)
    #sample_label = test_label[i]
    minIndex = np.argmin(disctance_vector)
    test_label[i] = train_label[minIndex]
    print ("The " + repr(i) +" sample, predicted label: " + repr(test_label[i]))
print (test_label)

发表评论

电子邮件地址不会被公开。 必填项已用*标注