# -*- coding: utf-8 -*-
'''
Created on 2018年1月18日
@author: Jason.F
@summary: 特征抽取-LDA方法,监督,发现最优化分类的特征子空间,基于特征呈正态分布和特征间相互独立
'''
import pandas as pd
import numpy as np
from sklearn.cross_validation import train_test_split
from sklearn.preprocessing import MinMaxScaler
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegression
import matplotlib.pyplot as plt
#第一步:导入数据,对原始d维数据集做标准化处理
df_wine = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/wine/wine.data',header=None)
df_wine.columns=['Class label','Alcohol','Malic acid','Ash','Alcalinity of ash','Magnesium','Total phenols','Flavanoids','Nonflavanoid phenols','Proanthocyanins','Color intensity','Hue','OD280/OD315 of diluted wines','Proline']
print ('class labels:',np.unique(df_wine['Class label']))
#分割训练集合测试集
X,y=df_wine.iloc[:,1:].values,df_wine.iloc[:,0].values
X_train,X_test,y_train,y_test=train_test_split(X,y,test_size=0.3,random_state=0)
#特征值缩放-标准化
stdsc=StandardScaler()
X_train_std=stdsc.fit_transform(X_train)
X_test_std=stdsc.fit_transform(X_test)
#第二步:对于每一类别计算d维的均值向量
np.set_printoptions(precision=4)
mean_vecs=[]
for label in range(1,4):
mean_vecs.append(np.mean(X_train_std[y_train==label],axis=0))
print ('MV %s: %s \n' %(label,mean_vecs[label-1]))
#第三步:构造类间的散布矩阵和类内的散布矩阵
d=13 #特征数量
#计算类内散布矩阵
#观察训练集的类别样本是否均匀,计算散布矩阵的前提是训练集的类标是均匀分布的
print ('class label distribution:%s' %np.bincount(y_train)[1:])
S_W=np.zeros((d,d))#初始化类内散布矩阵
for label,mv in zip(range(1,4),mean_vecs):
#class_scatter=np.zeros((d,d))
#for row in X[y==label]:
# row,mv =row.reshape(d,1),mv.reshape(d,1)
# class_scatter+= (row-mv).dot((row-mv).T)
#类标分布不均匀,对特征值做标准化,用标准化后的特征值计算散布矩阵
class_scatter=np.cov(X_train_std[y_train==label].T)#协方差矩阵是归一化的散布矩阵
S_W += class_scatter
print ('Within-class scatter matrix: %sx%s' %(S_W.shape[0],S_W.shape[1]))
#计算类间散布矩阵
mean_overall = np.mean(X_train_std,axis=0)
S_B=np.zeros((d,d))#初始化类间散布矩阵
for i ,mean_vec in enumerate(mean_vecs):
n=X_train_std[y_train==i+1,:].shape[0]
mean_vec=mean_vec.reshape(d,1)
mean_overall=mean_overall.reshape(d,1)
S_B+=n*(mean_vec-mean_overall).dot((mean_vec-mean_overall).T)
print ('Between-class scatter matrix: %sx%s' %(S_B.shape[0],S_B.shape[1]))
#第四部:计算类间类内乘积的矩阵的特征值和特征向量
eigen_vals,eigen_vecs=np.linalg.eig(np.linalg.inv(S_W).dot(S_B))
eigen_pairs=[(np.abs(eigen_vals[i]), eigen_vecs[:, i]) for i in range(len(eigen_vals))]
eigen_pairs=sorted(eigen_pairs,key=lambda k:k[0],reverse=True)
print ('Eigenvalues in decreasing order:\n')
for eigen_val in eigen_pairs:
print (eigen_val[0])
#可视化特征判定类别区分能力的图,按照特征值排序绘制出特征对线性判别信息保持程度
tot=sum(eigen_vals.real)
discr=[(i/tot) for i in sorted(eigen_vals.real,reverse=True)]
cum_discr=np.cumsum(discr)
plt.bar(range(1,14),discr,alpha=0.5,align='center',label='individual discriminability')
plt.step(range(1,14),cum_discr,where='mid',label='cumulative discriminability')
plt.ylabel('discriminability ratio')
plt.xlabel('Linear Discriminants')
plt.ylim([-0.1,1.1])
plt.legend(loc='best')
plt.show()
#第五步:选取前k个特征值所对应的特征向量,构造一个dXk维的转换矩阵W,其中特征向量以列的形式排列
w=np.hstack((eigen_pairs[0][1][:,np.newaxis].real,eigen_pairs[1][1][:,np.newaxis].real))#选取前2个特征,构建13X2维的映射矩阵W
print ('Matrix W:\n',w)
#第六步:使用转换矩阵W将样本映射到新的特征子空间
X_train_lda=X_train_std.dot(w)
X_test_lda=X_test_std.dot(w)
colors=['r','b','g']
markers=['s','x','o']
for l,c,m in zip(np.unique(y_train),colors,markers):
plt.scatter(X_train_lda[y_train == l, 0],X_train_lda[y_train == l, 1],c=c, label=l, marker=m)
plt.xlabel('LD 1')
plt.ylabel('LD 2')
plt.legend(loc='upper right')
plt.show()
#第五步:转换后的数据集进行线性训练
lr=LogisticRegression()
lr.fit(X_train_lda,y_train)
print ('Training accuracy:',lr.score(X_train_lda, y_train))
print ('Test accuracy:',lr.score(X_test_lda, y_test))# -*- coding: utf-8 -*-
'''
Created on 2018年1月1