开发: C++知识库 Java知识库 JavaScript Python PHP知识库人工智能区块链大数据移动开发嵌入式开发工具数据结构与算法开发测试游戏开发网络协议系统运维
教程: HTML教程 CSS教程 JavaScript教程 Go语言教程 JQuery教程 VUE教程 VUE3教程 Bootstrap教程 SQL数据库教程 C语言教程 C++教程 Java教程 Python教程 Python3教程 C#教程
数码: 电脑笔记本显卡显示器固态硬盘硬盘耳机手机 iphone vivo oppo 小米华为单反装机图拉丁

-> 人工智能 -> LDA(线性判别分析)详解 —— python -> 正文阅读

[人工智能]LDA(线性判别分析)详解 —— python

前言

正题

1.LDA的思想

?????????2. 瑞利商（Rayleigh quotient）与广义瑞利商（genralized Rayleigh quotient）?

前言

????????在主成分和因子分析中，我们对降维算法做了总结。这里我们就对另外一种经典的降维方法线性判别分析（Linear Discriminant Analysis, 以下简称LDA）做一个总结。LDA在模式识别领域（比如人脸识别，舰艇识别等图形图像识别领域）中有非常广泛的应用，因此我们有必要了解下它的算法原理。

　　　　在学习LDA之前，有必要将其自然语言处理领域的LDA区别开来，在自然语言处理领域， LDA是隐含狄利克雷分布（Latent Dirichlet Allocation，简称LDA），他是一种处理文档的主题模型。我们本文只讨论线性判别分析，因此后面所有的LDA均指线性判别分析。

在做具体解释之前，请允许我放上我之前的一些链接：

主成分分析 —— matlab ：传送门

主成分分析 —— python ：传送门

因子分析 —— matlab ：传送门

因子分析 —— python ：传送门

正题

1.LDA的思想

? ? ? ? 线性判别分析((Linear Discriminant Analysis?，简称 LDA)是一种经典的线性学习方法，在二分类问题上因为最早由 [Fisher,1936] 提出，亦称 ”Fisher 判别分析“。并且LDA也是一种监督学习的降维技术，也就是说它的数据集的每个样本都有类别输出。这点与主成分和因子分析不同，因为它们是不考虑样本类别的无监督降维技术。

????????LDA 的思想非常朴素：给定训练样例集，设法将样例投影到一条直线上，使得同样样例的投影尽可能接近、异样样例的投影点尽可能远离；在对新样本进行分类时，将其投影到同样的这条直线上，再根据投影点的位置来确定新样本的类别。其实可以用一句话概括：就是“投影后类内方差最小，类间方差最大”。

图为 LDA的二维示意图。”+“，”-“分别代表正侧和反侧，椭圆表示数据簇的外轮廓，虚线表示投影，红色实心圆和实心三角形分别表示两类样本投影后的中心点。

2. 瑞利商（Rayleigh quotient）与广义瑞利商（genralized Rayleigh quotient）?

? ? ? ? 我们先来看一下瑞丽商的定义。

? ? ? ? 瑞丽商是指这样的函数R(A,x)：

$R(A,x) = \frac{x^HAx}{x^Hx}$

? ? ? ? 其中x为非零向量，而A为 n*n 的Hermitan矩阵。所谓的Hermitan矩阵就是满足共轭转置矩阵和自己相等的矩阵，即 $A^{H}=A$ ?.?如果我们的矩阵A是实矩阵，则满足 $A^{T}=A$ 的矩阵即为Hermitan矩阵。

????????瑞利商R(A,x)有一个非常重要的性质，即它的最大值等于矩阵A最大的特征值，而最小值等于矩阵A的最小的特征值，也就是满足

$\lambda_{min} \leq \frac{x^HAx}{x^Hx} \leq \lambda_{max}$

? ? ? ? 至于证明过程就不在这里介绍了。当向量x是标准正交基时，即满足 $x^{H}x=1$ 时，瑞利商退化为： $R(A,x)=x^{H}Ax$ ，这个形式在谱聚类和PCA中都有出现。

以上就是瑞丽商的内容。

? ? ? ? 下面我们再来介绍一下广义的瑞丽商，

? ? ? ? 广义的瑞丽商是指这样的函数 R(A,B,x) ：

$R(A,x) = \frac{x^HAx}{x^HBx}$

????????其中x为非零向量，而A,B为n×n的Hermitan矩阵。B为正定矩阵。它的最大值和最小值是什么呢？其实我们只要通过将其通过标准化就可以转化为瑞利商的格式。我们令 $x=B^{-1/2}x'$ ,则分母转化为：

$x^HBx = x'^H(B^{-1/2})^HBB^{-1/2}x' = x'^HB^{-1/2}BB^{-1/2}x' = x'^Hx'$

而分子转化为：

$x^HAx = x'^HB^{-1/2}AB^{-1/2}x'$

此时我们的R(A,B,x)转化为R(A,B,x′)：

$R(A,B,x') = \frac{x'^HB^{-1/2}AB^{-1/2}x'}{x'^Hx'}$

????????利用前面的瑞利商的性质，我们可以很快的知道，R(A,B,x′)的最大值为矩阵 $B^{-1/2}AB^{-1/2}$ 的最大特征值，或者说矩阵 $B^{-1}A$ 的最大特征值，而最小值为矩阵 $B^{-1}A$ 的最小特征值。

3. 二类LDA原理

? ? ? ? 我们先介绍一下稍微简单，容易理解的二分类LDA入手，详细了解一下LDA的原理

? ? ? ? 首先给定数据集 D={(x1,y1),(x2,y2),...,((xm,ym))} ，其中任意样本xi为n维向量，yi∈{0,1}。另外我们定义Nj(j=0,1)为第j类样本的个数，Xj(j=0,1)为第j类样本的集合，而μj(j=0,1)为第j类样本的均值向量，定义Σj(j=0,1)为第j类样本的协方差矩阵（严格说是缺少分母部分的协方差矩阵）。

其中：

????????μj的表达式为：

$\mu_j = \frac{1}{N_j}\sum\limits_{x \in X_j}x\;\;(j=0,1)$

????????Σj的表达式为：

$\Sigma_j = \sum\limits_{x \in X_j}(x-\mu_j)(x-\mu_j)^T\;\;(j=0,1)$

我们将数据投影到一条直线上即可。我们假设我们的投影直线是向量w,则对任意一个样本本xi,它在直线w的投影为 $w^{T}x_{i}$ ,对于我们的两个类别的中心点μ0,μ1在在直线ww的投影为 $w^{T}u_{0}$ 和 $w^{T}u_{1}$ 。

由于LDA需要让不同类别的数据的类别中心之间的距离尽可能的大，也就是我们要最大化 $||w^T\mu_0-w^T\mu_1||_2^2$ ?,同时我们希望同一种类别数据的投影点尽可能的接近，也就是要同类样本投影点的协方差 $w^T\Sigma_0w$ ?和? $w^T\Sigma_1w$ ?尽可能的小，即最小化 $w^T\Sigma_0w+w^T\Sigma_1w$ ?。综上所述，我们的优化目标为：

$\underbrace{arg\;max}_w\;\;J(w) = \frac{||w^T\mu_0-w^T\mu_1||_2^2}{w^T\Sigma_0w+w^T\Sigma_1w} = \frac{w^T(\mu_0-\mu_1)(\mu_0-\mu_1)^Tw}{w^T(\Sigma_0+\Sigma_1)w}$

在这里，大家是否有很多问好？？？

就是 w ,? $w^T$ ?在哪，怎么算，下面就是我们求解的过程，在本小节最后就是哦！！！

我们一般定义类内散度矩阵Sw为：

$S_w = \Sigma_0 + \Sigma_1 = \sum\limits_{x \in X_0}(x-\mu_0)(x-\mu_0)^T + \sum\limits_{x \in X_1}(x-\mu_1)(x-\mu_1)^T$

同时定义类间散度矩阵Sb为：

$S_b = (\mu_0-\mu_1)(\mu_0-\mu_1)^T$

全局散度矩阵： $S_t=S_{b }+ S_{w}=\sum_{i=1}^{m} (x_i-\mu)(x_i-\mu)^T$

这样我们的优化目标重写为：

$\underbrace{arg\;max}_w\;\;J(w) = \frac{w^TS_bw}{w^TS_ww}$

????????仔细一看上式，这不就是我们的广义瑞利商嘛！这就简单了，利用我们第二节讲到的广义瑞利商的性质，我们知道我们的J(w')最大值为矩阵 $S^{-1/2}_wS_bS^{-1/2}_w$ 的最大特征值，而对应w'为 $S^{-1/2}_wS_bS^{-1/2}_w$ ?的最大特征值对应的特征向量! 而 $S{^{-1}_{w}}S_{b}$ 的特征值和 $S^{-1/2}_wS_bS^{-1/2}_w$ 的特征值相同,的特征向量 $S{^{-1}_{w}}S_{b}$ 的特征向量w和 $S^{-1/2}_wS_bS^{-1/2}_w$ 的特征向量w′满足 $w=S_{w}^{-1/2}w'$ 的关系!?

在这里我们就求到了 w 了哦！！！

? ? ? ? 注意到对于二类的时候， $S_bw$ 的方向恒平行于μ0?μ1,不妨令 $S_bw=\lambda (\mu_0-\mu _1)$ ，将其带入： $(S_w^{-1}S_b)w=\lambda w$ ，可以得到 $w=S_w^{-1}(\mu _0-\mu_1)$ ，也就是说我们只要求出原始二类样本的均值和方差就可以确定最佳的投影方向w了。

4.多类LDA原理

? ? ? ? 前面我们介绍了二分类的LDA，接下来我们再来看看多类别的LDA。

????????假设我们的数据集， $D=\{(x_1,y_1), (x_2,y_2), ...,((x_m,y_m))\}$ ，其中任意样本 $x_i$ ?为n维向量， $y_i \in \{C_1,C_2,...,C_k\}$ ?。我们定义 $N_j(j=1,2...k)$ 为j类样本的个数， $X_j(j=1,2...k)$ 为j类样本的集合，而 $\mu_j(j=1,2...k)$ ?为第j类样本的均值向量，定义 $\Sigma_j(j=1,2...k)$ ?为第j类样本的协方差矩阵。从在二类LDA里面定义的公式可以很容易的类推到多类LDA。

????????由于我们是多类向低维投影，则此时投影到的低维空间就不是一条直线，而是一个超平面了。假设我们投影到的低维空间的维度为d，对应的基向量为 $(w_1,w_2,...w_d)$ ，基向量组成的矩阵为W，它是一个n*d的矩阵。

????????此时我们的优化目标应该可以变成为:

$\frac{W^TS_bW}{W^TS_wW}$

其中， $S_b = \sum\limits_{j=1}^{k}N_j(\mu_j-\mu)(\mu_j-\mu)^T$ ， $\mu$ 为所有样本的均值向量。 $S_w = \sum\limits_{j=1}^{k}S_{wj} = \sum\limits_{j=1}^{k}\sum\limits_{x \in X_j}(x-\mu_j)(x-\mu_j)^T$

但是在这里会有一个问题？

????????就是 $W^TS_bW$ ?和 $W^TS_wW$ 都是矩阵，不是标量，无法作为一个标量函数来优化！也就是说，我们无法直接用二类LDA的优化方法，怎么办呢？

????????一般来说，我们可以用其他的一些替代优化目标来实现。

????????比如，常见的一个LDA多类优化目标函数定义为：

$\underbrace{arg\;max}_W\;\;J(W) = \frac{\prod\limits_{diag}W^TS_bW}{\prod\limits_{diag}W^TS_wW}$

? ? ? ? 其中， $\prod\limits_{diag}A$ ??为 A的主对角线元素的乘积，W为n×d的矩阵。

????????J(W)的优化过程可以转化为：

$J(W) = \frac{\prod\limits_{i=1}^dw_i^TS_bw_i}{\prod\limits_{i=1}^dw_i^TS_ww_i} = \prod\limits_{i=1}^d\frac{w_i^TS_bw_i}{w_i^TS_ww_i}$

此时，我再来观察一下上式，会发现最右边的式子，就是我们上面所讲的广义瑞丽商！！！最大值是矩阵 $S_w^{-1}S_b$ 的最大特征值，最大的d个值的乘积就是矩阵 $S_w^{-1}S_b$ ?的最大的d个特征值的乘积,此时对应的矩阵W为这最大的d个特征值对应的特征向量张成的矩阵。

　　　　由于W是一个利用了样本的类别得到的投影矩阵，因此它的降维到的维度d最大值为k-1。为什么最大维度不是类别数k呢？

????????因为Sb中每个μj?μ的秩为1，因此协方差矩阵相加后最大的秩为k(矩阵的秩小于等于各个相加矩阵的秩的和)，但是由于如果我们知道前k-1个μj后，最后一个μk可以由前k-1个μj线性表示，因此Sb的秩最大为k-1，即特征向量最多有k-1个。

5.LDA算法流程

输入：数据集 $D=\{(x_1,y_1), (x_2,y_2), ...,((x_m,y_m))\}$ ?，其中任意样本xixi为n维向量， $\large y_i \in \{C_1,C_2,...,C_k\}$ ，降维到的维度d。

输出：降维后的样本集D′

1) 计算类内散度矩阵Sw

2) 计算类间散度矩阵Sb?

3) 计算矩阵 $S^{-1}wSb$ ?

4）计算 $S^{-1}wSb$ 的最大的d个特征值和对应的d个特征向量(w1,w2,...wd)得到投影矩阵w

5)?对样本集中的每一个样本特征xi,转化为新的样本 $z_i=W^Tx_i$

6) 得到输出样本集? $D'=\{(z_1,y_1), (z_2,y_2), ...,((z_m,y_m))\}$ ?

二类LDApython举例：

1.首先看一下所需库

numpy 和 pandas

导入如下：

import numpy as np
import pandas as pd

2.读取数据

data = pd.read_excel('文件路径') #读取excel文件

data = pd.read_csv('文件路径') #读取csv文件

3.定义LDA函数

def LDA(x1, x2):
    n1 = np.shape(x1)[0]
    n2 = np.shape(x2)[0]
    m = n1
    print("m = ", m)
    x1_mean = np.mean(x1)
    x2_mean = np.mean(x2)
    print("x1_mean = ", x1_mean)
    print("x2_mean = ", x2_mean)
    sigma1 = np.dot((x1 - x1_mean),(x1 - x1_mean).T)
    sigma2 = np.dot((x2 - x2_mean),(x2 - x2_mean).T)
    Sw = sigma1 + sigma2
    print("Sw = ", Sw)
    w = Sw**(-1)*(x1_mean - x2_mean)
    return w