fisher线性判别算法python_降维算法-LDA线性判别分析实例

喜欢ヅ旅行 2022-11-06 02:59 229阅读 0赞

内容简介

线性判别分析LDA的基本概念

代码实例：第一部分使用python详细说明了LDA的计算过程； 第二部分记录了如何使用sklearn完成LDA。

什么是线性判别分析？

LDA，全名 Linear Discrimination Analysis， 是一种由监督学习的降维算法

LDA关心的是能够最大化类间的区分度的坐标轴成分。降特征投影到一个维度更小的k维子空间中，同时保持区分类别的信息。

原理：投影到维度更低的空间，是的投影后的点，会形成按类别区分一簇一簇的情况，相同类别的点，将会在投影后的空间中更接近。

用人话说，就是怎么样最简单地将数据分类。如左边的图，红色和蓝色在原来的坐标轴内不能很好地区分；因此重新画了一条坐标轴(带点的那条斜线)这条线上可以看到，红点和蓝点分别在两边；根据新的坐标轴，将数据投影到维度更低的空间中，即右图， 两种点分成了两个簇。

![4c66bf83559c2781d20907107eb752ef.png][]

课程中老师还举了这个游戏作为例子：我们的数据集就跟图中的铁坨坨一样，长得乱七八糟的，投影就是找到更适合分类的空间，是的数据集的信息可以被提取和理解；与PCA更关心方差，而LDA更关心分类。

![3a76d64b039b626f5caac67e6e3601ea.png][]

LDA计算过程理解

(这里涉及很多公式，简单复述一下自己的理解，如果有错或者描述不清楚的地方，请评论告诉我一声)

LDA目标：投影后的两类样本中心点尽量分离，即使类间距离最大化，

使用类间距离存在的问题： 如下图情况，按X1方向可以最大化两中心点的距离，但是分的不好，可见红圈和蓝圈在x1轴有很大的重叠部分；X2轴，虽然距离小，但是分类效果更好。

![4eabf50d6c7d2d132ee33f08b546ca97.png][]

散列值：样本点的密集程度，值越大，越分散，反之越集中； 即使得类内距离(类别里面个点的距离)尽量小。

python代码实例

依旧使用鸢尾花数据集来练习, 鸢尾花数据集导入方法有很多种；直接下载文档后用read\_csv();原来老师的代码是引用链接方法；个人觉得直接dataset 更方便，这里开头替换了一下取数的代码。在之前鸢尾花数据分析项目博客中有比较详细的解释这个方法，可以戳\[链接\]看看(https://blog.csdn.net/WHYbeHERE/article/details/106843722)

1.导入数据

\#导入模块

import numpy as np

import pandas as pd

from sklearn import datasets

iris = datasets.load\_iris() \#直接从datasets 导出数据

X, y = iris.data, iris.target \#切分数据集

df = pd.DataFrame(np.hstack((X, y.reshape(-1, 1))),index = range(X.shape\[0\]),columns=\['sepal\_length\_cm','sepal\_width\_cm','petal\_length\_cm','petal\_width\_cm','class'\] )

df.tail()

![269c1f754b7453f8748d86ff34ee2187.png][]

2.特征转换

from sklearn.preprocessing import LabelEncoder

\#这里原本是的X和y的赋值；被注释掉了，因为我们在上一步就分好了X和y

\#X = df\[\['sepal\_length\_cm','sepal\_width\_cm','petal\_length\_cm','petal\_width\_cm'\]\].values

\#y = df\['class'\].values

\#y原本的值是三种花的名称：'Setosa','Versicolor','Virginica'；需要将它们转换为电脑可以计算的数值。

enc = LabelEncoder()

label\_encoder = enc.fit(y)

y = label\_encoder.transform(y) + 1

\#label\_dict = \{1: 'Setosa', 2: 'Versicolor', 3:'Virginica'\}

![d7a382706ed89986156e0fc0a9f3ff5c.png][]

3. 求均值向量mi

接下来需要分别求三种鸢尾花数据在不同特征维度上的均值向量mi

![e3e466616fcd4dbfb5acd056052192d3.png][]

import numpy as np

np.set\_printoptions(precision=4)

mean\_vectors = \[\]

for cl in range(1,4):

mean\_vectors.append(np.mean(X\[y==cl\], axis=0))

print('Mean Vector class %s: %s\\n' %(cl, mean\_vectors\[cl-1\]))

结果输出：

![e6be782bcb3a0f64fc4a9d915a0ac9ac.png][]

每一种花都有四个值，因为我们有四个特征。

4.计算类内散步矩阵和类间散步矩阵

![7136e0765fa66eccac7e2fbfa415f46c.png][]

\#用代码实现上面的公式

S\_W = np.zeros((4,4))

for cl,mv in zip(range(1,4), mean\_vectors):

class\_sc\_mat = np.zeros((4,4)) \# scatter matrix for every class

for row in X\[y == cl\]:

row, mv = row.reshape(4,1), mv.reshape(4,1) \# make column vectors

class\_sc\_mat += (row-mv).dot((row-mv).T)

S\_W += class\_sc\_mat \# sum class scatter matrices

print('within-class Scatter Matrix:\\n', S\_W)

![e5509d54630d720e792813eefd273819.png][]

继续计算类间散布矩阵

![92daaf663eccee20b9b5fb5a85c2dcc7.png][]

overall\_mean = np.mean(X, axis=0)

S\_B = np.zeros((4,4))

for i,mean\_vec in enumerate(mean\_vectors):

n = X\[y==i+1,:\].shape\[0\]

mean\_vec = mean\_vec.reshape(4,1) \# make column vector

overall\_mean = overall\_mean.reshape(4,1) \# make column vector

S\_B += n \* (mean\_vec - overall\_mean).dot((mean\_vec - overall\_mean).T)

print('between-class Scatter Matrix:\\n', S\_B)

![6af9bd47a5e1558a2ebbd8f85220a8a8.png][]

5.求解矩阵的特征值：两个矩阵相乘

![904961fca9238d2af7ccb6f5f55ec122.png][]

eig\_vals, eig\_vecs = np.linalg.eig(np.linalg.inv(S\_W).dot(S\_B))

for i in range(len(eig\_vals)):

eigvec\_sc = eig\_vecs\[:,i\].reshape(4,1)

print('\\nEigenvector \{\}: \\n\{\}'.format(i+1, eigvec\_sc.real))

print('Eigenvalue \{:\}: \{:.2e\}'.format(i+1, eig\_vals\[i\].real))

结果得到四个特征值，对应四个特征向量

![f92de5e5e1405c1866ba5651d6b5c46e.png][]

6.特征值与特征向量：

特征向量：表示映射方向

特征值：特征向量的重要程度

\#Make a list of (eigenvalue, eigenvector) tuples

eig\_pairs = \[(np.abs(eig\_vals\[i\]), eig\_vecs\[:,i\]) for i in range(len(eig\_vals))\]

\# Sort the (eigenvalue, eigenvector) tuples from high to low

eig\_pairs = sorted(eig\_pairs, key=lambda k: k\[0\], reverse=True)

\# Visually confirm that the list is correctly sorted by decreasing eigenvalues

print('Eigenvalues in decreasing order:\\n')

for i in eig\_pairs:

print(i\[0\])

将特征值的结果排序：第一和第二个特征比较重要

![f53b378d79bb6a399caa7e55efc96c05.png][]

print('Variance explained:\\n')

eigv\_sum = sum(eig\_vals)

for i,j in enumerate(eig\_pairs):

print('eigenvalue \{0:\}: \{1:.2%\}'.format(i+1, (j\[0\]/eigv\_sum).real))

算一下百分比：结果也是选择前两个

![5ecc4df8a09540bfba9a20d3694b9abb.png][]

W = np.hstack((eig\_pairs\[0\]\[1\].reshape(4,1), eig\_pairs\[1\]\[1\].reshape(4,1)))

print('Matrix W:\\n', W.real)

两个特征向量，所以组合一起是一个4\*2的矩阵

![d805101a196d3d507c03ba66c30df107.png][]

X\_lda = X.dot(W)

assert X\_lda.shape == (150,2), "The matrix is not 150x2 dimensional."

\#画图

from matplotlib import pyplot as plt

def plot\_step\_lda():

ax = plt.subplot(111)

for label,marker,color in zip(

range(1,4),('^', 's', 'o'),('blue', 'red', 'green')):

plt.scatter(x=X\_lda\[:,0\].real\[y == label\],

y=X\_lda\[:,1\].real\[y == label\],

marker=marker,

color=color,

alpha=0.5,

label=label

)

plt.xlabel('LD1')

plt.ylabel('LD2')

leg = plt.legend(loc='upper right', fancybox=True)

leg.get\_frame().set\_alpha(0.5)

plt.title('LDA: Iris projection onto the first 2 linear discriminants')

\# hide axis ticks

plt.tick\_params(axis="both", which="both", bottom="off", top="off",

labelbottom="on", left="off", right="off", labelleft="on")

\# remove axis spines

ax.spines\["top"\].set\_visible(False)

ax.spines\["right"\].set\_visible(False)

ax.spines\["bottom"\].set\_visible(False)

ax.spines\["left"\].set\_visible(False)

plt.grid()

plt.tight\_layout

plt.show()

plot\_step\_lda()

两个坐标轴 LD 1和LD2 意味着我们将数据降到了二维

![b20040a5be7e4b69e4cee74bee4bca0b.png][]

\*直接使用sklearn中的LDA

from sklearn.discriminant\_analysis import LinearDiscriminantAnalysis as LDA

\# LDA

sklearn\_lda = LDA(n\_components=2)

X\_lda\_sklearn = sklearn\_lda.fit\_transform(X, y)

from matplotlib import pyplot as plt

def plot\_scikit\_lda(X, title):

ax = plt.subplot(111)

for label,marker,color in zip(

range(1,4),('^', 's', 'o'),('blue', 'red', 'green')):

plt.scatter(x=X\[:,0\]\[y == label\],

y=X\[:,1\]\[y == label\] \* -1, \# flip the figure

marker=marker,

color=color,

alpha=0.5,

label=label

)

plt.xlabel('LD1')

plt.ylabel('LD2')

leg = plt.legend(loc='upper right', fancybox=True)

leg.get\_frame().set\_alpha(0.5)

plt.title(title)

\# hide axis ticks

plt.tick\_params(axis="both", which="both", bottom="off", top="off",

labelbottom="on", left="off", right="off", labelleft="on")

\# remove axis spines

ax.spines\["top"\].set\_visible(False)

ax.spines\["right"\].set\_visible(False)

ax.spines\["bottom"\].set\_visible(False)

ax.spines\["left"\].set\_visible(False)

plt.grid()

plt.tight\_layout

plt.show()

plot\_step\_lda()

plot\_scikit\_lda(X\_lda\_sklearn, title='Default LDA via scikit-learn')

结果输出：虽然用sklearn的结果倒了过来，但是 对比两个图可以看出其实两个方法的结果都是一样的。

![87885ff67175f99cbbac052007cbf32d.png][]

![1acf97e4daa1ed05344a5e1922ed7aae.png][]

[4c66bf83559c2781d20907107eb752ef.png]: /images/20221023/28bd49664d1b45ebaa438f0ac783f391.png
[3a76d64b039b626f5caac67e6e3601ea.png]: /images/20221023/9715b42ee1cd46c5b27b44a585b3f82e.png
[4eabf50d6c7d2d132ee33f08b546ca97.png]: /images/20221023/a8eb92d4cab84e5fb7ab5c3c2fdc24f0.png
[269c1f754b7453f8748d86ff34ee2187.png]: /images/20221023/aef04adb467b4f41916fe389aff05505.png
[d7a382706ed89986156e0fc0a9f3ff5c.png]: /images/20221023/89acde6f9c0b450da08d17beea4b3517.png
[e3e466616fcd4dbfb5acd056052192d3.png]: /images/20221023/c7059cf465274701864ba53dd2c63f37.png
[e6be782bcb3a0f64fc4a9d915a0ac9ac.png]: /images/20221023/53849a773bf0477ba15ec8e66881f37c.png
[7136e0765fa66eccac7e2fbfa415f46c.png]: /images/20221023/796089df4fbb4d909dc9a7c8d2ad973d.png
[e5509d54630d720e792813eefd273819.png]: /images/20221023/1951a2d3c1c341b88b66d79a5f675350.png
[92daaf663eccee20b9b5fb5a85c2dcc7.png]: /images/20221023/3c9566b6a1264762aa858c34f1b71fbf.png
[6af9bd47a5e1558a2ebbd8f85220a8a8.png]: /images/20221023/615c4219bafe47b7ad76eec254387c63.png
[904961fca9238d2af7ccb6f5f55ec122.png]: /images/20221023/e4e99340b2bb465fa8f114d37483259b.png
[f92de5e5e1405c1866ba5651d6b5c46e.png]: /images/20221023/23af7599521d4051837b2aaf77ce9c8d.png
[f53b378d79bb6a399caa7e55efc96c05.png]: /images/20221023/4543989c0d3a4e778ac7caefa66ec274.png
[5ecc4df8a09540bfba9a20d3694b9abb.png]: /images/20221023/29d1fb92ef0d47c6bd6f0c311ef846c0.png
[d805101a196d3d507c03ba66c30df107.png]: /images/20221023/997b81aaa496417b8b7be092d4b8a581.png
[b20040a5be7e4b69e4cee74bee4bca0b.png]: /images/20221023/5624aa82836944b98e3b6f9976dc73bb.png
[87885ff67175f99cbbac052007cbf32d.png]: /images/20221023/8021302cc31142d8a9b83ebc68e47837.png
[1acf97e4daa1ed05344a5e1922ed7aae.png]: /images/20221023/3db96ea958614be99c942972d74cb36c.png