LDA(线性判别分析)及源码实现

最新推荐文章于 2024-02-12 21:32:02 发布

遗失的钥匙扣

最新推荐文章于 2024-02-12 21:32:02 发布

阅读量7.5k

点赞数

分类专栏：机器学习

机器学习专栏收录该内容

9 篇文章 0 订阅

订阅专栏

最近要用LDA对数据进行降维，下文是对网上写的较好的关于LDA文章的整理以便于日后复习。

原文链接：

http://www.cnblogs.com/zhangchaoyang/articles/2644095.html

http://blog.csdn.net/carson2005/article/details/8652586

1、判别分析

首先搞清楚什么叫判别分析？Discriminant Analysis就是根据研究对象的各种特征值判别其类型归属问题的一种多变量统计分析方法。根据判别标准不同，可以分为距离判别、Fisher判别、Bayes判别法等。比如在KNN中用的就是距离判别，朴素贝叶斯分类用的就是Bayes判别法。本文要讲的线性判别分析就是用是Fisher判别式。根据判别函数的形式，可以分为线性判别和非线性判别。

2、什么是LDA

线性判别式分析（Linear Discriminant Analysis），简称为LDA。也称为Fisher线性判别（Fisher Linear Discriminant，FLD），是模式识别的经典算法，在1996年由Belhumeur引入模式识别和人工智能领域。

基本思想是将高维的模式样本投影到最佳鉴别矢量空间，以达到抽取分类信息和压缩特征空间维数的效果，投影后保证模式样本在新的子空间有最大的类间距离和最小的类内距离，即模式在该空间中有最佳的可分离性。

3、LDA与PCA

LDA与PCA都是为了在对原始数据降维之后进行分类。PCA主要是从特征的协方差角度，去找到比较好的投影方式。LDA更多的是考虑了标注，即希望投影后不同类别之间数据点的距离更大，同一类别的数据点更紧凑。PCA是无监督的方式，它没有分类标签，降维之后需要采用K-Means或自组织映射网络等无监督的算法进行分类。LDA是有监督的方式，它先对训练数据进行降维，然后找出一个线性判别函数。

下图是一个LDA和PCA区别的例子：

两个类别，class1的点都是圆圈，class2的点都是十字。图中有两条直线，斜率在1左右的这条直线是PCA选择的映射直线，斜率在 -1左右的这条直线是LDA选择的映射直线。其余不在这两条直线上的点是原始数据点。可以看到由于LDA考虑了“类别”这个信息（即标注），映射后，可以很好的将class1和class2的点区分开。

4、LDA计算过程

4.1 两类线性判别分析

给定N个特征为d维的样例 $x^{(i)}\{x_1^{(i)},x_2^{(i)},...,x_d^{(i)}\}$ ，其中有N₁个样例属于类别w₁，另外N₂个样例属于类别w₂。现在我们要将原始数据降低到只有一维，降维函数（或者叫投影函数）是： $y=w^Tx$ ，最后我们就依靠每个样例对应的y值来判别它属于哪一类。

形象的图求如下

我们就是要找到这个最佳的w，使得样例映射到y后最易于区分。

定义每类样例的均值点： $u_i=\frac{1}{N_i}\sum_{x\in{w_i}}x$

样例投影到y后有均值点为： $\tilde{u_i}=\frac{1}{N_i}\sum_{x\in{w_i}}{w^Tx}=w^Tu_i$

我们希望投影后两类样例中心尽量地分离，即

$|\tilde{u_1}-\tilde{u_2}|=|w^T(u_1-u_2)|$ 越大越好。

同时我们希望投影之后类内部的方差 $\tilde{s_i}^2=\sum_{y\in{w_i}}(y-\tilde{u_i})^2$ 越小越好。

由于得到我们的目标函数：

$max\qquad{J(w)}=\frac{|\tilde{u_1}-\tilde{u_2}|^2}{\tilde{s_1}^2+\tilde{s_2}^2}\qquad(1)$

又是个最优化问题。最终解得

$w=(s_1+s_2)^{-1}(u_1-u_2)$ ，s₁和s₂分别中原始样例的方差。

这就是Fisher在1936年提出的线性判别分析。

如果 $y=w^Tx-w^Tu>0$ （u是所有样本的均值），就属于类别C₁，否则就属于类别C₂。

实际上使用线性回归得到的直线方向就是二值分类中LDA求得的直线方向。

4.2 多类线性判别分析

假设有C个类别，降以一维已经不能满足分类要求了，我们需要k个基向量来做投影，W=[w₁|w₂|...|w_k] 。样本点在这k维投影后的结果为[y₁,y₂,...,y_k]，且有 $y_i=w_i^Tx\qquad{y}=W^Tx$ 。

$S_w=\sum_{i=1}^C{S_{wi}}$

$S_{wi}=\sum_{x\in{w_i}}(x-u_i)(x-u_i)^T$

$S_B=\sum_{i=1}^C{N_i(u_i-u)(u_i-u)^T}$

$u=\frac{1}{N}\sum_{\forall{x}}x=\frac{1}{N}\sum_{x\in{w_i}}N_iu_i$

同样是求一个类似于(1)式的最优化问题，我们得到

$S_w^{-1}S_Bw_i=\lambda{w_i}$

即w_i是矩阵 $S_w^{-1}S_B$ 的特征向量。首先求出 $S_w^{-1}S_B$ 的特征值，然后取其前k个特征向量组成W即可，因为特征值大的对应的特征向量分割性好。另外S_B的秩至多为C-1，所以 $S_w^{-1}S_B$ 不为0的特征值至多有C-1个，所以k最大为C-1。

由于 $S_w^{-1}S_B$ 不一定是对称矩阵，因此k个特征向量不一定正交，这也是与PCA不同的地方。求 $S_w^{-1}S_B$ 的特征向量不能采用奇异值分解的方式，而因该采用更通用的求一般方阵特征向量的方式。特征值的求法有很多，求一个D * D的矩阵的时间复杂度是O(D^3), 也有一些求Top k的方法，比如说幂法，它的时间复杂度是O(D^2 * k)。

那降维之后又如何根据y值来判别分类呢？取[y₁,y₂,...,y_k]中最大的那个就是所属的分类。这么说对于有C个类别的分类问题，我们最多只能分出C-1个类别来？

5、LDA的一些说明

第一，降维后的维度是多少？

PCA降维是直接和数据维度相关的，比如原始数据是n维的，那么PCA后，可以任意选取1维、2维，一直到n维都行（当然是对应特征值大的那些）。LDA 降维是直接和类别的个数相关的，与数据本身的维度没关系，比如原始数据是n维的，一共有C个类别，那么LDA降维之后，一般就是1维，2维到C-1维进行选择（当然对应的特征值也是最大的一些），举个例子，假设图象分类，两个类别正例反例，每个图象10000维特征，那么LDA之后，就只有1维特征，并且这维特征的分类能力最好。
PS：对于很多两类分类的情况，LDA之后就剩下1维，找到分类效果最好的一个阈值貌似就可以了。

第二，投影的坐标系是否正交

PCA投影的坐标系都是正交的，而LDA根据类别的标注，关注分类能力，因此不保证投影到的坐标系是正交的（一般都不正交）。

6、自己写的LDA程序

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Drawing;
using Emgu.CV;
using Emgu.CV.Structure;
using Emgu.CV.UI;
using Emgu.Util;

namespace DecreaseFeatureLDA
{
    public class LDA
    {
        Matrix<float> sample;//样本矩阵
        int[] label;//样本标签
        int classNum;//类别数
       
        
       public LDA(float[,] s, int[] l,int num)
        {
            sample = new Matrix<float>(s);
            label = l;
            classNum = num;
           
        }
       public void DoLDA(out Matrix<float> eivects,out Matrix<float>eivalues)
        {
            int sampleNum = sample.Rows;//样本数
            int sampleDim = sample.Cols;//维数
            
            List<int> labelIndex = new List<int>(classNum);//记录每个类的样本的索引
            for ( int i = 0; i < label.Length; )
            {
                
                int temp = 0;
                for ( int j = i; j < label.Length; j++)
                {
                   
                    if (label[i] == label[j])
                    {
                       
                        temp++;
                        if (j==label.Length-1)
                        {
                            labelIndex.Add(temp);
                            i = label.Length;
                            break;
                        }
                        
                    }
                    else
                    {
                        labelIndex.Add(temp);
                        i = j;
                        break;
                    }

                }
            }
            labelIndex.Insert(0, 0);
           
            //计算类中心
            Matrix<float> classCenter = new Matrix<float>(classNum, sampleDim);
            int start=0;
            //求每个类的均值
            for (int i = 0; i < classNum; i++)
            {
                start += labelIndex.ElementAt(i);
                for (int j = 0; j < labelIndex.ElementAt(i+1); j++)
                {
                    for (int k1 = 0; k1 < sampleDim; k1++)
                    {
                        classCenter[i, k1] += sample[j+start, k1];
                    }
                }
                for (int k2 = 0; k2 < sampleDim; k2++)
                {
                    classCenter[i, k2] /= labelIndex.ElementAt(i + 1);
                }
               
            }
            
        
          //计算各类样本的协方差矩阵Sw
            Matrix<float>[] covSw=new Matrix<float>[classNum];
            Matrix<float> Sw = new Matrix<float>(sampleDim, sampleDim);
            start = 0;
            for (int i = 0; i < classNum; i++)
            {
                covSw[i] = new Matrix<float>(sampleDim, sampleDim);
                start += labelIndex.ElementAt(i);
                for (int j = 0; j < labelIndex.ElementAt(i + 1); j++)
                {
                    for (int k1 = 0; k1 < sampleDim;k1++)
                    {
                        for (int k2 = 0; k2 < sampleDim; k2++)
                        {
                            covSw[i][k1, k2] += (sample[j + start, k1] - classCenter[i, k1]) * (sample[j + start, k2] - classCenter[i, k2]);
                        }

                    }
                }
                Sw += covSw[i];
            }
            
            //计算各类间协方差矩阵
            Matrix<float> sampleCenter = new Matrix<float>(1, sampleDim);//样本中心
        
            for (int i = 0; i < sampleNum; i++)
            {
                for (int j = 0; j < sampleDim; j++)
                {
                    sampleCenter[0, j] += sample[i, j];
                }
            }
            for (int j = 0; j < sampleDim; j++)
            {
                sampleCenter[0, j] /= sampleNum;
            }
            //求类间协方差Sb
            Matrix<float>[] covSb = new Matrix<float>[classNum];
            Matrix<float> Sb=new Matrix<float> (sampleDim,sampleDim);
            for (int i = 0; i < classNum; i++)
            {
                covSb[i] = new Matrix<float>(sampleDim, sampleDim);
                for (int k1 = 0; k1 < sampleDim; k1++)
                {
                    for (int k2 = 0; k2 < sampleDim; k2++)
                    {
                        covSb[i][k1, k2] = labelIndex.ElementAt(i+1)*(classCenter[i, k1] - sampleCenter[0, k1]) * (classCenter[i, k2] - sampleCenter[0, k2]);
                    }
                }
            }
            for (int i = 0; i < classNum; i++)
            {
                Sb += covSb[i]; 
            }
            Matrix<float> W=new Matrix<float> (sampleDim,sampleDim);
            Matrix<float> Sw_inverse=new Matrix<float> (sampleDim,sampleDim);
            Matrix<float> V=new Matrix<float> (sampleDim,sampleDim);
            CvInvoke.cvInvert(Sw.Ptr,Sw_inverse,Emgu.CV.CvEnum.SOLVE_METHOD.CV_LU);
             Matrix<float> S=new Matrix<float> (sampleDim,sampleDim);
             S = Sw_inverse * Sb;
            eivects=new Matrix<float> (sampleDim,sampleDim);
            eivalues=new Matrix<float> (sampleDim,sampleDim);
            Matrix<float> U=new Matrix<float> (sampleDim,sampleDim);
            CvInvoke.cvSVD(S.Ptr, eivalues.Ptr, U.Ptr, eivects.Ptr, Emgu.CV.CvEnum.SVD_TYPE.CV_SVD_V_T);
           
        }
    }
}

如训练数据为{{ 2.95,6.63 }, { 2.53,7.79 }, { 3.57,5.65}, { 3.16,5.47}, {2.58,4.46 }, { 2.16,6.22 } ,{3.27,3.52}} 时，对应的类别标签为{ 1,1, 1, 1,0,0,0 }，代入程序：

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using Emgu.CV;
using Emgu.CV.Structure;
using Emgu.CV.UI;
using Emgu.Util;

namespace DecreaseFeatureLDA
{
    class Program
    {
        static void Main(string[] args)
        {
            float[,] source = new float[7, 2] { { 2.95f,6.63f }, { 2.53f,7.79f }, { 3.57f,5.65f}, { 3.16f,5.47f}, {2.58f,4.46f }, { 2.16f,6.22f } ,{3.27f,3.52f}};
            int[] label = new int[7] { 1,1, 1, 1,0,0,0 };
            float r = 0.9f;
            Matrix<float> eivects;
            Matrix<float> eivalues;
            LDA lda = new LDA(source, label,2);
            lda.DoLDA(out eivects,out eivalues);
            Matrix<float> e = new Matrix<float>(2, 1);
            e[0, 0] = eivects[0, 1];
            e[1, 0] = eivects[1, 1];
            Matrix<float> S=new Matrix<float> (source);
            Matrix<float> project = S * e;                   
        }
    }
}

特征值eivalues=[18,147,0],选择特征为18.147对应的特征向量[0.974,0.226]进行投影，投影后project=[4.370,4.222,4.753,4.313,3.520,3.508,3.980],可见在投影后的轴上点4处能将两类数据完全分开。

遗失的钥匙扣

关注

0
点赞
踩
18

收藏

觉得还不错? 一键收藏
0
评论
LDA(线性判别分析)及源码实现

最近要用LDA对数据进行降维，下文是对网上写的较好的关于LDA文章的整理以便于日后复习。1、LDA首先搞清楚判别分析？Discriminant Analysis就是根据研究对象的各种特征值判别其类型归属问题的一种多变量统计分析方法。线性判别式分析（Linear Discriminant Analysis），简称为LDA。也称为Fisher线性判别（Fisher L
复制链接

扫一扫