LDA与PCA算法与应用

来源:互联网 发布:哪个挂号软件好用抢号 编辑:程序博客网 时间:2024/04/30 10:05

【原文:http://blog.csdn.net/whuheming/article/details/19809889】

(一)LDA算法

LDA的全称是Linear Discriminant Analysis(线性判别分析),一种supervised learning有些资料上也称为是Fisher’s Linear Discriminant,因为它被Ronald Fisher发明自1936年,Discriminant这次词我个人的理解是,一个模型,不需要去通过概率的方法来训练、预测数据,比如说各种贝叶斯方法,就需要获取数据的先验、后验概率等等。

LDA的原理是,将带上标签的数据(点),通过投影的方法,投影到维度更低的空间中,使得投影后的点,会形成按类别区分,一簇一簇的情况,相同类别的点,将会在投影后的空间中更接近。线性判别分析的基本思想是将高维的模式样本投影到最佳鉴别矢量空间,即把高维空间中的数据点投影到一条直线去,将多维降为一维,并且要求投影后各样本的类间散布距离最大,同时类内散布距离最小。即投影之后的数据点更具有线性分类性

LDA是一种线性分类器。对于K-分类的一个分类问题,会有K个线性函数:

201101081455464493.jpg

当满足条件:对于所有的j,都有Yk > Yj,的时候,我们就说x属于类别k。对于每一个分类,都有一个公式去算一个分值,在所有的公式得到的分值中,找一个最大的,就是所属的分类了。

假设用来区分二分类的直线(投影函数)为:

 

类别i的原始中心点为:

201101081455478264.jpg(Di表示属于类别i的点)

    类别i投影后的中心点为:
201101081455488231.jpg
   衡量类别i投影后,类别点之间的分散程度(方差)为:
201101081455487326.jpg
   最终我们可以得到一个下面的公式,表示LDA投影到w后的损失函数,也称为准则函数

201101081455484785.jpg

我们分类的目标是找到一个最优化的W,使得类别内的点距离越近越好(集中),类别间的点越远越好。上式准则函数分母表示每一个类别内的方差之和,方差越大表示一个类别内的点越分散,分子为两个类别各自的中心点的距离的平方,我们最大化J(w)就可以求出最优的w。想要求出最优的w,可以使用拉格朗日乘子法,但是现在我们得到的J(w)里面,w是不能被单独提出来的,我们就得想办法将w单独提出来。

我们定义样本类内离散度矩阵和总类内离散度矩阵如下所示:

201101081455498656.jpg

代入上式准则函数,可以简化为:

201101081455505982.jpg  这样就可以用最喜欢的拉格朗日乘子法了,但是还有一个问题,如果分子、分母是都可以取任意值的,那就会使得有无穷解,我们将分母限制为长度为1(这是用拉格朗日乘子法一个很重要的技巧,在下面将说的PCA里面也会用到,如果忘记了,请复习一下高数),并作为拉格朗日乘子法的限制条件,代入得到:

201101081455516963.jpg 

这同样是一个求特征值的问题,我们求出的第i大的特征向量,就是对应的Wi了。

(二) PCA算法

主成分分析法(PCA)是多元统计分析中用来分析数据的一种方法,它是用一种较少数量的特征对样本进行描述以达到降低特征空间维数的方法,它的本质实际上是K-L变换。PCA是一种unsupervised learning,最著名的应用应该是在人脸识别中特征提取及数据维,我们知道输入200*200大小的人脸图像,单单提取它的灰度值作为原始特征,则这个原始特征将达到40000维,这给后面分类器的处理将带来极大的难度。著名的人脸识别Eigenface算法就是采用PCA算法,用一个低维子空间描述人脸图像,同时用保存了识别所需要的信息。

主成分分析(PCA)的原理就是将一个高维向量x,通过一个特殊的特征向量矩阵U,投影到一个低维的向量空间中,表征为一个低维向量y,并且仅仅损失了一些次要信息。也就是说,通过低维表征的向量和特征向量矩阵,可以基本重构出所对应的原始高维向量。在人脸识别中,特征向量矩阵U称为特征脸(eigenface)空间,因此其中的特征向量ui进行量化后可以看出人脸轮廓,在下面的实验中可以看出。以人脸识别为例,说明下PCA的应用。假设有N个人脸训练样本,每个样本由其像素灰度值组成一个向量xi,则样本图像的像素点数即为xi的维数,M=width*height ,由向量构成的训练样本集为:{x1,x2,.....xN}

该样本集的平均向量(又称为平均脸)为:

样本集的协方差矩阵为:

求出协方差矩阵的特征向量ui和对应的特征值,这些特征向量组成的矩阵U就是人脸空间的正交基底,用它们的线性组合可以重构出样本中任意的人脸图像,并且图像信息集中在特征值大的特征向量中,即使丢弃特征值小的向量也不会影响图像质量。
将协方差矩阵的特征值按大到小排序:。由大于对应的特征向量构成主成分,主成分构成的变换矩阵为:
这样每一幅人脸图像都可以投影到构成的特征脸子空间中,U的维数为M×d。有了这样一个降维的子空间,任何一幅人脸图像都可以向其作投影,即并获得一组坐标系数,即低维向量y,维数d×1,为称为KL分解系数。这组系数表明了图像在子空间的位置,从而可以作为人脸识别的依据。
以下是我采用openCV写的基于PCA人脸识别的代码:
[html] view plaincopy在CODE上查看代码片派生到我的代码片
  1. //***************************** Face Recognize ***********************************  
  2. void train(){     
  3.     // load training data    
  4.     nTrainFaces = loadFaceImgArray("faceSample\\train.txt");    
  5.     if( nTrainFaces < 2){    
  6.         fprintf(    
  7.             stderr,    
  8.             "Need 2 or more training faces\n"    
  9.             "Input file contains only %d\n",    
  10.             nTrainFaces       
  11.         );    
  12.         return;    
  13.     }    
  14.     // do PCA on the training faces    
  15.     doPCA();    
  16.     // project the training images onto the PCA subspace    
  17.     projectedTrainFaceMat = cvCreateMat(nTrainFaces, nEigens, CV_32FC1);    
  18.     for(int i = 0; i < nTrainFaces; i ++){    
  19.         cvEigenDecomposite(    
  20.             faceImgArr[i],    
  21.             nEigens,    
  22.             eigenVectArr,    
  23.             0, 0,    
  24.             pAvgTrainImg,    
  25.             projectedTrainFaceMat->data.fl + i*nEigens    
  26.         );    
  27.     }    
  28.     // store the recognition data as an xml file    
  29.     storeTrainingData();    
  30. }   
  31.   
  32. int loadFaceImgArray(char* filename){    
  33.     FILE* imgListFile = 0;    
  34.     char imgFilename[512];    
  35.     int iFace, nFaces = 0;    
  36.   
  37.     // open the input file    
  38.     imgListFile = fopen(filename, "r");    
  39.         
  40.     // count the number of faces    
  41.     while( fgets(imgFilename, 512, imgListFile) ) ++ nFaces;    
  42.     rewind(imgListFile);    
  43.     
  44.     // allocate the face-image array and person number matrix    
  45.     faceImgArr = (IplImage **)cvAlloc( nFaces*sizeof(IplImage *) );    
  46.     personNumTruthMat = cvCreateMat( 1, nFaces, CV_32SC1 );    
  47.     
  48.     // store the face images in an array    
  49.     for(iFace=0; iFace<nFaces; iFace++){    
  50.     //read person number and name of image file    
  51.     fscanf(imgListFile, "%d %s", personNumTruthMat->data.i+iFace, imgFilename);    
  52.     
  53.     // load the face image    
  54.     faceImgArr[iFace] = cvLoadImage(imgFilename, CV_LOAD_IMAGE_GRAYSCALE);    
  55.     }    
  56.       
  57.     fclose(imgListFile);    
  58.     
  59.     return nFaces;    
  60. }    
  61.   
  62. void doPCA(){    
  63.     int i;    
  64.     CvTermCriteria calcLimit;    
  65.     CvSize faceImgSize;    
  66.     
  67.     // set the number of eigenvalues to use    
  68.     nEigens = nTrainFaces - 1;    
  69.     
  70.     // allocate the eigenvector images    
  71.     faceImgSize.width = faceImgArr[0]->width;    
  72.     faceImgSize.height = faceImgArr[0]->height;    
  73.     eigenVectArr = (IplImage**)cvAlloc(sizeof(IplImage*) * nEigens);    
  74.     for(i=0; i<nEigens; i++){    
  75.     eigenVectArr[i] = cvCreateImage(faceImgSize, IPL_DEPTH_32F, 1);    
  76.     }    
  77.     // allocate the eigenvalue array    
  78.     eigenValMat = cvCreateMat( 1, nEigens, CV_32FC1 );    
  79.     // allocate the averaged image    
  80.     pAvgTrainImg = cvCreateImage(faceImgSize, IPL_DEPTH_32F, 1);      
  81.     // set the PCA termination criterion    
  82.     calcLimit = cvTermCriteria( CV_TERMCRIT_ITER, nEigens, 1);     
  83.     // compute average image, eigenvalues, and eigenvectors    
  84.     cvCalcEigenObjects(    
  85.     nTrainFaces,    
  86.     (void*)faceImgArr,    
  87.     (void*)eigenVectArr,    
  88.     CV_EIGOBJ_NO_CALLBACK,    
  89.     0,    
  90.     0,    
  91.     &calcLimit,    
  92.     pAvgTrainImg,    
  93.     eigenValMat->data.fl    
  94.     );    
  95. }  
  96.   
  97. void storeTrainingData(){    
  98.     CvFileStorage* fileStorage;    
  99.     int i;    
  100.        
  101.     // create a file-storage interface    
  102.     fileStorage = cvOpenFileStorage( "facedata.xml", 0, CV_STORAGE_WRITE);    
  103.     
  104.     // store all the data    
  105.     cvWriteInt( fileStorage, "nEigens", nEigens);    
  106.     cvWriteInt( fileStorage, "nTrainFaces", nTrainFaces );    
  107.     cvWrite(fileStorage, "trainPersonNumMat", personNumTruthMat, cvAttrList(0, 0));    
  108.     cvWrite(fileStorage, "eigenValMat", eigenValMat, cvAttrList(0,0));    
  109.     cvWrite(fileStorage, "projectedTrainFaceMat", projectedTrainFaceMat, cvAttrList(0,0));    
  110.     cvWrite(fileStorage, "avgTrainImg", pAvgTrainImg, cvAttrList(0,0));    
  111.     
  112.     for(i=0; i<nEigens; i++){    
  113.     char varname[200];    
  114.     sprintf( varname, "eigenVect_%d", i);    
  115.     cvWrite(fileStorage, varname, eigenVectArr[i], cvAttrList(0,0));    
  116.     }    
  117.     
  118.     //release the file-storage interface    
  119.     cvReleaseFileStorage( &fileStorage );    
  120. }    
  121. // Face Recognize by PCA EigenFace...  
  122. int PCA_recognize(IplImage* faceRegion)  
  123. {  
  124.     CvMat* trainPersonNumMat = 0;   // the person numbers during training        
  125.     float* projectedTestFace = 0;  
  126.   
  127.      // load the saved training data    
  128.     if( !loadTrainingData( &trainPersonNumMat ) ) return -1;  
  129.   
  130.      // project the test face region onto the PCA subspace    
  131.     projectedTestFace = (float*)cvAlloc( nEigens*sizeof(float) );  
  132.     cvEigenDecomposite(    
  133.         faceRegion,    
  134.         nEigens,    
  135.         eigenVectArr,    
  136.             0, 0,    
  137.         pAvgTrainImg,    
  138.         projectedTestFace    
  139.     );  
  140.     int iNearest, nearest;   
  141.     iNearest = findNearestNeighbor(projectedTestFace);     
  142.     nearest = trainPersonNumMat->data.i[iNearest];   
  143.     return nearest;  
  144. }  
0 0
原创粉丝点击