线性判别分析（LDA）原理详解及实现

线性判别分析（LDA）的实现步骤

线性判别分析（Linear Discriminant Analysis，LDA）是一种用于分类和降维的方法。其目标是寻找一个线性组合，使得不同类别的数据在新的维度上具有最大的可分性。

LDA的基本步骤

1. 计算类内散度矩阵和类间散度矩阵：

   • 类内散度矩阵 $S_W$ 表示同一类别中数据点的分散程度。
   • 类间散度矩阵 $S_B$ 表示不同类别之间的分散程度。

2. 计算矩阵 $S_W^{-1}{S}_B$ 的特征值和特征向量：

   • 这一步的目的是找到一个矩阵，当数据点乘以这个矩阵后，最大化类间散度同时最小化类内散度。

3. 选择特征值最大的特征向量：

   • 选取最大的几个特征值对应的特征向量。这些特征向量决定了最佳的新坐标轴。

4. 将数据投影到新的坐标轴上：

   • 使用选定的特征向量将原始数据投影到新的坐标系，实现降维。

数学公式

类内散度矩阵：

$$S_W=\sum _{i=1}^c\sum _{\mathbf{x}\in D_i}(\mathbf{x}-{\mathbf{m}}_i)(\mathbf{x}-{\mathbf{m}}_i{)}^T$$

其中，$c$ 是类别数，$D_i$ 是第 $i$ 个类别的数据集，${\mathbf{m}}_i$ 是第 $i$ 个类别的均值向量。

类内散度矩阵 $S_W$ 描述了同一类别内部数据点的分散程度。以下是一个具体的例子来说明这一点。

示例

假设我们有两个类别的数据集，每个类别包含若干二维数据点。

• 类别 1 的数据点为：$\{(1,2),(2,3),(3,3)\}$，均值向量为 ${\mathbf{m}}_1=(2,2.67)$
• 类别 2 的数据点为：$\{(6,8),(7,8),(8,9)\}$，均值向量为 ${\mathbf{m}}_2=(7,8.33)$

计算类内散度矩阵

类内散度矩阵 $S_W$ 由两个类别的散度矩阵 $S_1$ 和 $S_2$ 相加得到：

$$S_W=S_1+S_2$$

其中，$S_1$ 和 $S_2$ 分别是类别 1 和类别 2 的散度矩阵：

$$S_1=\sum _{\mathbf{x}\in \text{类别1}}(\mathbf{x}-{\mathbf{m}}_1)(\mathbf{x}-{\mathbf{m}}_1{)}^T$$

$$S_2=\sum _{\mathbf{x}\in \text{类别2}}(\mathbf{x}-{\mathbf{m}}_2)(\mathbf{x}-{\mathbf{m}}_2{)}^T$$

例如，对于类别 1 中的一个点 $(1,2)$，进行外积运算，其贡献到 $S_1$ 的部分是：

$$\begin{array}{rl}(\mathbf{x}-{\mathbf{m}}_1)(\mathbf{x}-{\mathbf{m}}_1{)}^T& =\left[\begin{array}{c}1-2\\ 2-2.67\end{array}\right]\left[\begin{array}{cc}1-2& 2-2.67\end{array}\right]\\ & =\left[\begin{array}{c}-1\\ -0.67\end{array}\right]\left[\begin{array}{cc}-1& -0.67\end{array}\right]\\ & =\left[\begin{array}{cc}1& 0.67\\ 0.67& 0.45\end{array}\right]\end{array}$$

类似地，我们会计算类别 1 中其他点的贡献，以及类别 2 中所有点的贡献，然后将这些矩阵相加得到 $S_W$。

为什么用外积

外积运算能生成一个方阵，这个方阵的每个元素代表了数据点在各个维度上的分散情况。具体来说：

• 方阵对角线上的元素表示数据点在对应维度上的分散程度（方差）。
• 非对角线元素表示不同维度间的协方差，即这些维度间分散程度的相关性。

因此，通过外积运算，我们能够获得一个全面的描述数据点分散情况的矩阵。

代码

import numpy as np

# 假设数据集
class1 = np.array([[1, 2], [2, 3], [3, 3]])  # 类别1的数据点
class2 = np.array([[6, 8], [7, 8], [8, 9]])  # 类别2的数据点

# 计算每个类别的均值向量
mean1 = np.mean(class1, axis=0)
mean2 = np.mean(class2, axis=0)

# 计算类内散度矩阵
S_W = np.zeros((2, 2))  # 初始化类内散度矩阵
for x in class1:
    S_W += np.outer((x - mean1), (x - mean1))
for x in class2:
    S_W += np.outer((x - mean2), (x - mean2))

print("类内散度矩阵 S_W:", S_W)

类间散度矩阵：

$$S_B=\sum _{i=1}^c{N}_i({\mathbf{m}}_i-\mathbf{m})({\mathbf{m}}_i-\mathbf{m}{)}^T$$

其中，$N_i$ 是第 $i$ 个类别中的样本数量，$\mathbf{m}$ 是所有样本的全局均值向量。

类间散度矩阵（$S_B$）用于描述不同类别之间的数据点的分散程度。它关注的是类别的均值向量之间的差异，而不是类别内部的数据点。

示例

假设我们有两个类别的数据集，每个类别包含若干二维数据点。

• 类别 1 的数据点为：$\{(1,2),(2,3),(3,3)\}$，均值向量为 ${\mathbf{m}}_1=(2,2.67)$
• 类别 2 的数据点为：$\{(6,8),(7,8),(8,9)\}$，均值向量为 ${\mathbf{m}}_2=(7,8.33)$

计算类间散度矩阵

类间散度矩阵计算的是不同类别的均值向量之间的差异。对于上面的例子，我们首先计算总体均值向量 $\mathbf{m}$：

$$\mathbf{m}=\text{所有数据点的均值}$$

然后，对于每个类别，我们计算它的均值向量与总体均值向量之间的差，然后进行外积运算，并乘以该类别中的样本数量。例如，对于类别 1：

$$S_{B1}=N_1({\mathbf{m}}_1-\mathbf{m})({\mathbf{m}}_1-\mathbf{m}{)}^T$$

类似地，我们计算类别 2 的贡献 $S_{B2}$。类间散度矩阵 $S_B$ 是这些矩阵的总和：

$$S_B=S_{B1}+S_{B2}$$

为什么用外积

使用外积计算类间散度矩阵可以帮助我们捕捉不同类别间均值向量的差异。外积运算产生的矩阵描述了这些差异在各个维度上的分布情况。

为什么用样本加权

从数学的角度来看，这种加权方式保证了在优化过程中，大类别的特征不会被小类别的特征所掩盖。
在计算特征值和特征向量时，这种加权方法有助于突出那些在数据集中占据主导地位的类别特征，从而在降维后的特征空间中实现有效的类别分离。

代码

import numpy as np

# 假设数据集
class1 = np.array([[1, 2], [2, 3], [3, 3]])  # 类别1的数据点
class2 = np.array([[6, 8], [7, 8], [8, 9]])  # 类别2的数据点

# 计算每个类别的均值向量和总体均值向量
mean1 = np.mean(class1, axis=0)
mean2 = np.mean(class2, axis=0)
overall_mean = np.mean(np.concatenate((class1, class2)), axis=0)

# 计算类间散度矩阵
S_B = np.zeros((2, 2))
S_B += len(class1) * np.outer((mean1 - overall_mean), (mean1 - overall_mean))
S_B += len(class2) * np.outer((mean2 - overall_mean), (mean2 - overall_mean))

print("类间散度矩阵 S_B:", S_B)

矩阵 $S_W^{-1}{S}_B$

在线性判别分析（LDA）中，计算 $S_W^{-1}{S}_B$ 的特征值和特征向量是核心步骤。这一步骤的目的是为了找到最佳的数据投影方向，以最大化类间散度同时最小化类内散度。

$S_W^{-1}{S}_B$ 的含义

• $S_W^{-1}$：这是类内散度矩阵 $S_W$ 的逆矩阵。它的作用是“调整”类间散度矩阵 $S_B$，在数学上，逆矩阵可以被视为一种“反效应”或“补偿”机制。$S_W^{-1}$ 实际上是在减小那些在原始空间中类内分散较大的特征的权重。
• $S_B$：类间散度矩阵，它描述了不同类别之间的分散程度。
• $S_W^{-1}{S}_B$：将 $S_B$ 与 $S_W^{-1}$ 相乘，其效果是在保持类间分离的同时，调整（或权衡）类内的分散程度。

特征值和特征向量

• 特征值：$S_W^{-1}{S}_B$ 的特征值表示了在对应特征向量定义的方向上，数据投影后类间散度与类内散度的比率。
• 特征向量：这些特征向量定义了数据从原始高维空间到新低维空间的投影方向。

为什么要求特征值和特征向量

• 在LDA中，我们的目标是寻找一个可以最大化类间分离度的新空间。这可以通过选择那些使得 $S_W^{-1}{S}_B$ 的特征值最大化的特征向量来实现。
• 特征值越大，表示在对应的特征向量方向上，类别间的分散程度相对于类别内的分散程度越大。因此，这些方向是区分不同类别最有效的。
• 换句话说，这些特征向量定义了一个新的坐标系，在这个坐标系中，不同类别的数据最大限度地分离开来。

通过计算 $S_W^{-1}{S}_B$ 的特征值和特征向量，我们能够确定数据降维的最佳方向，这有助于提高分类的准确性和效率。

代码

import numpy as np

# 假设 S_W 是类内散度矩阵，S_B 是类间散度矩阵

# 计算 S_W 的逆矩阵
S_W_inv = np.linalg.inv(S_W)

# 计算乘积矩阵 S_W^{-1}S_B
S_W_inv_S_B = np.dot(S_W_inv, S_B)

# 计算特征值和特征向量
eigenvalues, eigenvectors = np.linalg.eig(S_W_inv_S_B)

# 对特征值进行排序，并选择对应的特征向量
indices = np.argsort(eigenvalues)[::-1]
top_eigenvectors = eigenvectors[:, indices[:n_components]]

LDA（线性判别分析）算法的Python实现详解

线性判别分析（LDA）是一种常用于降维和分类的统计方法。以下是对LDA算法的Python实现进行详细解释。

算法步骤和代码解释

1. 计算类别均值和总体均值

• 对于数据集中的每个类别，计算该类别所有数据点的均值向量。
• 同时，计算整个数据集的总体均值向量。

2. 计算类内散布矩阵 $S_W$

• 类内散布矩阵 $S_W$ 量化了同一类别内数据点的分散程度。
• 对于每个类别，计算数据点与该类别均值向量的差，进行外积运算，并累加所有类别的结果。

3. 计算类间散布矩阵 $S_B$

• 类间散布矩阵 $S_B$ 描述了不同类别之间的分散程度。
• 对每个类别，计算其均值向量与总体均值向量的差，进行外积运算，并乘以该类别中的样本数量，再累加所有类别的结果。

4. 求解特征值和特征向量

• 使用矩阵 $S_W^{-1}{S}_B$ 的特征值和特征向量来确定最佳的投影方向。
• 使用伪逆（np.linalg.pinv(S_W)）来处理 $S_W$ 可能不可逆的情况。
• 求解 $S_W^{-1}{S}_B$ 的特征值和特征向量。

5. 选择主成分

• 根据特征值的大小，选择最重要的 num_components 个特征向量作为投影方向。

6. 正交化鉴别向量

• 如果需要正交化，使用QR分解（np.linalg.qr）对特征向量进行正交化处理。
• 正交化后的鉴别向量用于降维或分类。

代码

#X代表数据特征数组，y代表类别（标签），num_components代表降维数量
def lda(X, y, num_components,orth=1):
    class_means = {}
    for label in np.unique(y):
        class_means[label] = np.mean(X[y == label], axis=0)
    overall_mean = np.mean(X, axis=0)

    # 类内散布矩阵 S_W
    S_W = np.zeros((X.shape[1], X.shape[1]))
    for label, mean_vec in class_means.items():
        S_W += np.dot((X[y == label] - mean_vec).T, (X[y == label] - mean_vec))

    # 类间散布矩阵 S_B
    S_B = np.zeros((X.shape[1], X.shape[1]))
    for label, mean_vec in class_means.items():
        n = X[y == label].shape[0]
        mean_diff = (mean_vec - overall_mean).reshape(X.shape[1], 1)
        S_B += n * (mean_diff).dot(mean_diff.T)

    eigenvalues, eigenvectors = np.linalg.eigh(np.linalg.pinv(S_W).dot(S_B))
    idx = np.argsort(eigenvalues)[::-1]
    eigenvectors = eigenvectors[:, idx]
    eigenvectors = eigenvectors[:, :num_components]
    # 正交化鉴别向量
    if orth==0:
        print("unorth vector")
        return eigenvectors.T
    else :
        q, r = np.linalg.qr(eigenvectors)
        return q.T  # 返回正交化的鉴别向量

部分代码详解

class_means[label] = np.mean(X[y == label], axis=0)

这行代码是线性判别分析（LDA）算法中的一个关键步骤，它用于计算每个类别的均值向量。具体来说，这行代码的作用和意义如下：

计算类别均值向量

• 代码功能：对于数据集中的每个类别，计算该类别所有数据点的均值向量。
• 实现方式：使用 np.mean 函数沿着特定的轴（axis=0）计算均值。

数学原理

• 对于一个给定的类别（标记为 label），这行代码首先筛选出所有属于这个类别的数据点（X[y == label]），然后计算这些点在每个特征维度上的平均值。

• 数学表示为：

  $${\mathbf{m}}_{\text{label}}=\frac{1}{N_{\text{label}}}\sum _{\mathbf{x}\in D_{\text{label}}}\mathbf{x}$$

  其中，$D_{\text{label}}$ 是属于类别 label 的数据点集合，$N_{\text{label}}$ 是这个集合中的数据点数量，${\mathbf{m}}_{\text{label}}$ 是计算得到的均值向量。

S_W += np.dot((X[y == label] - mean_vec).T, (X[y == label] - mean_vec))

这段代码是线性判别分析（LDA）算法中用于计算类内散度矩阵 $S_W$ 的关键部分。以下是对这部分代码的详细解释。

代码功能

这部分代码的主要目的是计算类内散度矩阵 $S_W$，它度量了数据集中每个类别内部数据点的分散程度。

实现步骤

1. 迭代每个类别

• 代码首先通过 for label, mean_vec in class_means.items(): 循环遍历每个类别。
• label 是类别标签，mean_vec 是该类别的均值向量。

2. 计算差异矩阵

• 对于每个类别，计算所有数据点与该类别均值向量的差异，即 X[y == label] - mean_vec。
• 这里 X[y == label] 筛选出属于当前类别 label 的所有数据点，mean_vec 是该类别的均值向量。

3. 外积运算和累加

• 接着，计算差异向量的外积 np.dot((X[y == label] - mean_vec).T, (X[y == label] - mean_vec))。
• 外积运算结果是一个矩阵，描述了当前类别内部的数据点相对于均值的分散程度。
• 然后将这个结果累加到 S_W 中，以便对所有类别执行相同的操作。

数学原理

• 类内散度矩阵的计算公式是：

  $$S_W=\sum _{i=1}^c\sum _{\mathbf{x}\in D_i}(\mathbf{x}-{\mathbf{m}}_i)(\mathbf{x}-{\mathbf{m}}_i{)}^T$$

  其中，$c$ 是类别数，$D_i$ 是第 $i$ 个类别的数据集，${\mathbf{m}}_i$ 是第 $i$ 个类别的均值向量。

• 这个计算过程实际上是在评估每个类别内部数据点相对于其均值向量的分散情况。