Python 中的 Fisher 线性判别式

Question

我有 Fisher 线性判别式，我需要用它来将作为高维矩阵的示例 A 和 B 减少到简单的 2D，这与 LDA 完全一样，每个示例都有 A 和 B 类，因此如果我有第三个例如，它们也有 A 类和 B 类，第四、第五和 n 个例子总是有 A 类和 B 类，因此我想简单地使用 Fisher 线性判别式将它们分开。我对机器学习很陌生，所以我不知道如何分开我的课程，我一直在通过眼睛和编码来遵循公式。从我正在阅读的内容来看，我需要对我的数据应用线性变换，以便为它找到一个好的阈值，但首先我需要找到最大化函数。对于这样的任务，我设法找到了 Sw 和 Sb，但我不知道如何从那里开始......

我还需要找到最大化函数的地方。

该最大化函数给了我一个特征值解决方案：

我为每个类所拥有的是 2 个示例的 5x2 矩阵。例如：

Example 1 
Class_A = [
201, 103,
40, 43,
23, 50,
12, 123,
99, 78
]
Class_B = [
   201, 129,
   114, 195,
   180, 90,
   69, 62,
   76, 90
]

Example 2   
Class_A = [
68, 98,
201, 203,
78, 212,
49, 5,
204, 78
]   
Class_B = [
   52, 19,
   220, 219,
   159, 195,
   99, 23,
   46, 50
]

我尝试为上面的示例找到 Sw，如下所示：

Example_1_Class_A = np.dot(Example_1_Class_A,  np.transpose(Example_1_Class_A))
Example_1_Class_B = np.dot(Example_1_Class_B,  np.transpose(Example_1_Class_B))

Example_2_Class_A = np.dot(Example_2_Class_A,  np.transpose(Example_2_Class_A))
Example_2_Class_B = np.dot(Example_2_Class_B,  np.transpose(Example_2_Class_B))

Sw = sum([Example_1_Class_A, Example_1_Class_B, Example_2_Class_A, Example_2_Class_B], axis=0)

至于 Sb，我试过这样：

Example_1_Class_A_mean = Example_1_Class_A.mean(axis=0)
Example_1_Class_B_mean = Example_1_Class_B.mean(axis=0)
         
Example_2_Class_A_mean = Example_2_Class_A.mean(axis=0)
Example_2_Class_B_mean = Example_2_Class_B.mean(axis=0)
         
Example_1_Class_A_Sb = np.dot(Example_1_Class_A_mean, np.transpose(Example_1_Class_A_mean))
Example_1_Class_B_Sb = np.dot(Example_1_Class_B_mean, np.transpose(Example_1_Class_B_mean))
         
Example_2_Class_A_Sb = np.dot(Example_2_Class_A_mean, np.transpose(Example_2_Class_A_mean))
Example_2_Class_B_Sb = np.dot(Example_2_Class_B_mean, np.transpose(Example_2_Class_B_mean))
Sb = sum([Example_1_Class_A_Sb, Example_1_Class_B_Sb, Example_2_Class_A_Sb, Example_2_Class_B_Sb], axis=0) 

问题是，我不知道我的 Sw 和 Sb 还能做什么，我完全迷失了。基本上，我需要做的是从这里到这里：

对于给定的示例 A 和示例 B，我是否仅针对类 As 和仅针对类 b 分离集群

Answer 1

在回答你的问题之前，我先来了解一下 PCA 和 (F)LDA 的基本区别。在 PCA 中，您对基础类一无所知，但您假设有关类可分离性的信息在于数据的方差。因此，您旋转原始轴（有时称为将所有数据投影到新轴上），这样您的第一个新轴指向方差最大的方向，第二个新轴垂直于第一个并指向大多数剩余方差，等等。这样，PCA 变换会产生与原始空间具有相同维度的（子）空间。比你只能取前 2 个维度，拒绝其余的维度，因此从k 中得到降维 尺寸仅为 2。

LDA 的工作方式略有不同。在这种情况下，您事先知道数据中有多少个类，并且可以找到它们的均值和协方差矩阵。什么Fisher准则确实发现在哪些类之间的平均被最大化，而一个方向在同一时间总变异性最小化（总变异性是一个意思类内协方差矩阵的）。每两个类只有一个这样的行。这就是为什么当你的数据有C类时，LDA 最多可以为你提供C-1维度，而不管原始数据的维度。在您的情况下，这意味着由于您只有 2 个类 A 和 B，您将获得一维投影，即一条线。这正是您图片中的内容：原始二维数据投影到一条线上。线的方向是本征问题的解。让我们生成与您的图片类似的数据：

a = np.random.multivariate_normal((1.5, 3), [[0.5, 0], [0, .05]], 30)
b = np.random.multivariate_normal((4, 1.5), [[0.5, 0], [0, .05]], 30)
plt.plot(a[:,0], a[:,1], 'b.', b[:,0], b[:,1], 'r.')
mu_a, mu_b = a.mean(axis=0).reshape(-1,1), b.mean(axis=0).reshape(-1,1)
Sw = np.cov(a.T) + np.cov(b.T)
inv_S = np.linalg.inv(Sw)
res = inv_S.dot(mu_a-mu_b)  # the trick
####
# more general solution
#
# Sb = (mu_a-mu_b)*((mu_a-mu_b).T)
# eig_vals, eig_vecs = np.linalg.eig(inv_S.dot(Sb))
# res = sorted(zip(eig_vals, eig_vecs), reverse=True)[0][1] # take only eigenvec corresponding to largest (and the only one) eigenvalue
# res = res / np.linalg.norm(res)

plt.plot([-res[0], res[0]], [-res[1], res[1]]) # this is the solution
plt.plot(mu_a[0], mu_a[1], 'cx')
plt.plot(mu_b[0], mu_b[1], 'yx')
plt.gca().axis('square')

# let's project data point on it
r = res.reshape(2,)
n2 = np.linalg.norm(r)**2
for pt in a:
    prj = r * r.dot(pt) / n2
    plt.plot([prj[0], pt[0]], [prj[1], pt[1]], 'b.:', alpha=0.2)
for pt in b:
    prj = r * r.dot(pt) / n2
    plt.plot([prj[0], pt[0]], [prj[1], pt[1]], 'r.:', alpha=0.2)

使用针对两类问题的巧妙技巧计算得到的投影。您可以在第 1.6 节中阅读有关它的详细信息。

关于您在问题中提到的“示例”。我相信您需要为每个示例重复该过程，因为它是一组不同的数据点，可能具有不同的分布。此外，请注意估计的均值 (mu_a, mu_b) 和类协方差矩阵与生成数据的矩阵略有不同，尤其是对于小样本量。

Answer 2

数学

有关更多信息，请参阅https://sebastianraschka.com/Articles/2014_python_lda.html#lda-in-5-steps。

使用 Iris 实现

由于您想使用 LDA 进行降维但仅提供二维数据，因此我将展示如何在 iris 数据集上执行此过程。

让我们导入库

    import pandas as pd
    import numpy as np
    import sklearn as sk
    from collections import Counter
    from sklearn import datasets
    
    # load dataset and transform to pandas df
    X, y = datasets.load_iris(return_X_y=True)
    X = pd.DataFrame(X, columns=[f'feat_{i}' for i in range(4)])
    y = pd.DataFrame(y, columns=['labels'])
    tot = pd.concat([X,y], axis=1)
    # calculate class means
    class_means = tot.groupby('labels').mean()
    total_mean = X.mean()

由class_means以下给出：

class_means

    feat_0  feat_1  feat_2  feat_3
labels              
0   5.006   3.428   1.462   0.246
1   5.936   2.770   4.260   1.326
2   6.588   2.974   5.552   2.026

为此，我们首先从每个观察中减去类均值（基本上我们x - m_i从上面的等式计算）。从每个观察中减去相应的类平均值。因为我们要计算

x_mi = tot.transform(lambda x: x - class_means.loc[x['labels']], axis=1).drop('labels', 1)

def kronecker_and_sum(df, weights):
    S = np.zeros((df.shape[1], df.shape[1]))
    for idx, row in df.iterrows():
        x_m = row.as_matrix().reshape(df.shape[1],1)
        S += weights[idx]*np.dot(x_m, x_m.T)
    return S

# Each x_mi is weighted with 1. Now we use the kronecker_and_sum function to calculate the within-class scatter matrix S_w
S_w = kronecker_and_sum(x_mi, 150*[1])

mi_m = class_means.transform(lambda x: x - total_mean, axis=1)
# Each mi_m is weighted with the number of observations per class which is 50 for each class in this example. We use kronecker_and_sum to calculate the between-class scatter matrix.

S_b=kronecker_and_sum(mi_m, 3*[50])

eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(S_w).dot(S_b))

我们只需要考虑与零显着不同的特征值（在这种情况下只有前两个）

eig_vals
array([ 3.21919292e+01,  2.85391043e-01,  6.53468167e-15, -2.24877550e-15])

X用对应于最高特征值的两个特征向量的矩阵进行变换

eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(S_w).dot(S_b))

我们将维度从 4 减少到 2，并以这样的方式选择空间，使类可以很好地分开。

Scikit-learn 用法

Scikit 也有 LDA 支持。我们在几十行中所做的事情可以通过以下几行代码来完成：

from sklearn import discriminant_analysis
lda = discriminant_analysis.LinearDiscriminantAnalysis(n_components=2)
X_trafo_sk = lda.fit_transform(X,y)
pd.DataFrame(np.hstack((X_trafo_sk, y))).plot.scatter(x=0, y=1, c=2, colormap='viridis')

我没有在这里给出图表，因为它与我们的派生示例中的相同（除了 180 度旋转）。