算法原理

紧接着上文：https://blog.csdn.net/just_sort/article/details/102690086 。在上文我们使用LDA实现了一个二分类任务。那么数据有大于2种类别，假设为C类，这时候怎么办呢？在上文我们定义的“类间散度矩阵”
$S_{b} = (μ_{1} - μ_{2}) (μ_{1} - μ_{2})^{T}$
就不再适用，所以我们这里引入“全局散度矩阵”：
$S_{t} = S_{w} + S_{b} = \sum_{i = 1}^{n} (x_{i} - μ) (x_{i} - μ)^{T}$
其中 $μ$ 是所有示例的均值向量。将类内散度矩阵 $S_{w}$ 重定义为每个类别的散度矩阵之和，即：
$S_{w} = \sum_{i = 1}^{N} S_{w_{i}}$ ，其中：
$S_{w_{i}} = \sum_{x \in X_{i}} (x - μ_{i}) (x - μ_{i})^{T}$ 。从上面两个式子推得:
$S_{b} = S_{t} - S_{w} = \sum_{i = 1}^{N} m_{i} (μ_{i} - μ) (μ_{i} - μ)^{T}$ ，其中N代表类别数， $m_{i}$ 代表第i类的示例数。
关于这个等式的推导可以看南瓜书，我这里截图过来一下：

这里偷了一张图，可以更好的理解这个算法。
显然，多分类LDA有多种实现方法：使用 $S_{b}, S_{w}, S_{t}$ 三者中的任意两个即可。常见一种实现是采用优化目标：
$\max_{W} \frac{t r (W^{T} S_{b} W)}{t r (W^{T} S_{w} W)}$
式3.35就是我们上篇博客写的“广义瑞利商”。其中的tr()为矩阵的迹，一个n×n的对角矩阵A的主对角线（从左上方至右下方的对角线）上各个元素的总和被称为矩阵A的迹（或迹数），一般记作tr(A)。
这个优化目标实际上等价于求解多个w组合成W，那么该问题就等价于求解多个上一章的优化目标，使用相同的方法，可以求得下式：
$S_{b} W = λ S_{w} W$
即是： $S_{w}^{- 1} S_{b} W = λ W$
W的闭式解为 $S_{w}^{- 1} S_{b}$ 的 $d^{'}$ 个最大非零广义特征值对应的特征向量组成的矩阵, $d^{'} < = N - 1$ 。
如果将 $W$ 视为一个投影矩阵，则多分类LDA将样本投影到 $d^{'}$ 维空间， $d^{'}$ 通常远小于数据的原有属性 $d^{'}$ 。于是，可以通过这个投影来减小样本点的维数，且投影过程中使用了类别信息，因此LDA也常常被视为一种经典的监督降维技术。

接下来使用sklearn中的aris数据，并使用LDA算法对其进行降维，并可视化。

代码

#coding=utf-8
import numpy as np
from sklearn.datasets import load_iris
import matplotlib.pyplot as plt
# 这是sklearn中实现的LDA，待会我们会比较自己实现的LDA和它的区别
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis

# k为目标
def LDA(X, y, k):
    label_ = list(set(y))
    X_classify = {}
    for label in label_:
        X1 = np.array([X[i] for i in range(len(X)) if y[i] == label])
        X_classify[label] = X1

    miu = np.mean(X, axis=0)
    miu_classify = {}
    for label in label_:
        miu1 = np.mean(X_classify[label], axis=0)
        miu_classify[label] = miu1

    # St = np.dot((X - mju).T, X - mju)
    # 计算类内散度矩阵Sw
    Sw = np.zeros((len(miu), len(miu)))
    for i in label_:
        Sw += np.dot((X_classify[i] - miu_classify[i]).T, X_classify[i] - miu_classify[i])

    #Sb = St-Sw
    # 计算类内散度矩阵Sb
    Sb = np.zeros((len(miu), len(miu)))
    for i in label_:
        Sb += len(X_classify[i]) * np.dot((miu_classify[i] - miu).reshape(
            (len(miu), 1)), (miu_classify[i] - miu).reshape((1, len(miu))))

    # 计算S_w^{-1}S_b的特征值和特征矩阵
    eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(Sw).dot(Sb))
    sorted_indices = np.argsort(eig_vals)
    # 提取前k个特征向量
    topk_eig_vecs = eig_vecs[:, sorted_indices[:-k - 1:-1]]
    return topk_eig_vecs

def main():
    iris = load_iris()
    X = iris.data
    y = iris.target

    W = LDA(X, y, 2)
    X_new = np.dot(X, W)
    plt.scatter(X_new[:, 0], X_new[:, 1], marker='o', c=y)
    plt.show()

    # 和sklearn的函数对比
    lda = LinearDiscriminantAnalysis(n_components=2)
    lda.fit(X, y)
    X_new = lda.transform(X)
    plt.scatter(X_new[:, 0], X_new[:, 1], marker='o', c=y)
    plt.show()


main()

降到二维的效果图

(我们实现的LDA算法的降维结果)
可以看到使用LDA算法成功实现了多分类数据的降维。

参考资料

《周志华机器学习》
https://blog.csdn.net/z962013489/article/details/79918758

机器学习算法之线性判别分析(LDA多分类)

算法原理

代码

降到二维的效果图

参考资料