Accord Vision: Linear Discriminant Analysis in C#

Linear discriminant analysis (LDA) is a method used in statistics and machine learning to find a linear combination of features which best characterize or separates two or more classes of objects or events. The resulting combination may be used as a linear classifier, or, more commonly, for dimensionality reduction before later classification.

Download source code

Download sample applications

Browse the source code online.

The code presented here is also part of the Accord.NET Framework. The Accord.NET Framework is a framework for developing machine learning, computer vision, computer audition, statistics and math applications. It is based on the already excellent AForge.NET Framework. Please see the starting guide for mode details. The latest version of the framework includes the latest version of this code plus many other statistics and machine learning tools.

Motivations

The goals of LDA are somewhat similar to those of PCA. But different from LDA, PCA is an unsupervised technique and as such does not include label information of the data, effectively ignoring this often useful information. For instance, consider two cigar-like clusters in 2 dimensions, one cigar having y = c and the other y = –c (with c being a arbitrary constant), as the image below suggests:

This example was adapted from the
note on Fisher Linear Discriminant Analysis by Max Welling.

The cigars are positioned in parallel and very closely together, such that the variance in the total dataset, ignoring the labels, is in the direction of the cigars. For classification, this would be a terrible projection, because all labels will get mixed and we will destroy any useful information:

A much more useful projection is orthogonal to the cigars, which would perfectly separate the data-cases:

The first row of images was obtained by performing PCA on the original data. The left image is the PCA projection using the first two components. However, as the analysis says the first component accounts for 96% of the information, one can infer all other components could be safely discarded. The result is the image on the right, clearly a not very useful projection.

The second row of images was obtained by performing LDA on the original data. The left image is the LDA projection using the first two dimensions. However, as the analysis says the first dimension accounts for 100% of the information, all other dimensions could be safely discarded. Well, this is actually true, as we can see on the result in the right image.

Analysis Overview

Linear Discriminant Analysis is closely related to principal component analysis (PCA) and factor analysis in that both look for linear combinations of variables which best explain the data. LDA explicitly attempts to model the difference between the classes of data. PCA on the other hand does not take into account any difference in class, and factor analysis builds the feature combinations based on differences rather than similarities. Linear Discriminant Analysis considers maximizing the following objective:

$J(w) = \frac{w^T S_B w}{w^T S_W w}$

where

$S_B = \sum_{c=1}^C (\mu_c - \bar{x}) (\mu_c - \bar{x})^T$

$S_W = \sum_{c=1}^C \sum_{i \in c} (x_i - \mu_c) (x_i - \mu_c)^T$

are the Between-Class Scatter Matrix and Within-Class Scatter Matrix, respectively. The optimal solution can be found by computing the Eigenvalues of S_BS_W^-1 and taking the Eigenvectors corresponding to the largest Eigen values to form a new basis for the data. Those can be easily obtained by computing the generalized eigenvalue decomposition of S_B and S_W.

Source Code

The source code below shows the main algorithm for Linear Discriminant Analysis.

/// <summary>
///   Computes the Multi-Class Linear Discriminant Analysis algorithm.
/// </summary>
public virtual void Compute()
{
    // Compute entire data set measures
    Means = Tools.Mean(source);
    StandardDeviations = Tools.StandardDeviation(source, totalMeans);
    double total = dimension;

    // Initialize the scatter matrices
    this.Sw = new double[dimension, dimension];
    this.Sb = new double[dimension, dimension];


    // For each class
    for (int c = 0; c < Classes.Count; c++)
    {
        // Get the class subset
        double[,] subset = Classes[c].Subset;
        int count = subset.GetLength(0);

        // Get the class mean and number of samples
        double[] mean = Tools.Mean(subset);


        // Continue constructing the Within-Class Scatter Matrix
        double[,] Swi = Tools.Scatter(subset, mean);
        Sw = Sw.Add(Swi);

        // Continue constructing the Between-Class Scatter Matrix
        double[] d = mean.Subtract(totalMeans);
        double[,] Sbi = d.Multiply(d.Transpose()).Multiply(total/count);
        Sb = Sb.Add(Sbi);

        // Store some additional information
        this.classScatter[c] = Swi;
        this.classCount[c] = count;
        this.classMeans[c] = mean;
        this.classStdDevs[c] = Tools.StandardDeviation(subset, mean);
    }



    // Compute eigen value decomposition
    EigenValueDecomposition evd = new EigenValueDecomposition(Matrix.Inverse(Sw).Multiply(Sb));

    // Gets the eigenvalues and corresponding eigenvectors
    double[] evals = evd.RealEigenValues;
    double[,] eigs = evd.EigenVectors;


    // Sort eigen values and vectors in ascending order
    eigs = Matrix.Sort(evals, eigs, new GeneralComparer(ComparerDirection.Descending, true));


    // Calculate translation bias
    this.bias = eigs.Transpose().Multiply(totalMeans).Multiply(-1.0);

    // Store information
    this.EigenValues = evals;
    this.DiscriminantMatrix = eigs;


    // Create projections
    this.result = new double[dimension, dimension];
    for (int i = 0; i < dimension; i++)
        for (int j = 0; j < dimension; j++)
            for (int k = 0; k < dimension; k++)
                result[i, j] += source[i, k] * eigenVectors[k, j];


    // Computes additional information about the analysis and creates the
    //  object-oriented structure to hold the discriminants found.
    createDiscriminants();
}

/// <summary>Projects a given matrix into discriminant space.</summary>
/// <param name="matrix">The matrix to be projected.</param>
/// <param name="dimensions">The number of dimensions to consider.</param>
public virtual double[,] Transform(double[,] data, int dimensions)
{
    int rows = data.GetLength(0);
    int cols = data.GetLength(1);

    double[,] r = new double[rows, dimensions];

    // multiply the data matrix by the selected eigenvectors
    for (int i = 0; i < rows; i++)
        for (int j = 0; j < dimensions; j++)
            for (int k = 0; k < cols; k++)
                r[i, j] += data[i, k] * eigenVectors[k, j];

    return r;
}

Using The Code

Code usage is simple, as it follows the same object model in the previous PCA example.

   // Creates the Linear Discriminant Analysis of the given source
   lda = new LinearDiscriminantAnalysis(data, labels);


   // Compute the Linear Discriminant Analysis
   lda.Compute();


   // Perform the transformation of the data using two dimensions
   double[,] result = lda.Transform(data, 2);

Considerations

Besides being useful in many situations, in many practical cases linear discriminants are not suitable. A nice insight is that LDA can be extended for use in non-linear classification via the kernel trick. Using kernels, the original observations can be effectively mapped into a higher dimensional non-linear space. Linear classification in this non-linear space is then equivalent to non-linear classification in the original space.

While the code available here already contains code for Kernel Discriminant Analysis, this is something I'll address in the next post. If you have any suggestions or questions, please leave me a comment.

Finally, as usual, I hope someone finds this information useful.

References:

Welling, Max; Fisher Linear Discriminant Analysis.

Availabe in: http://www.ics.uci.edu/~welling/classnotes/papers_class/Fisher-LDA.pdf.

Gutierrez-Osuna, Ricardo; Introduction to Pattern Analysis.

Available in: http://research.cs.tamu.edu/prism/lectures/pr/pr_l10.pdf.

Accord Vision

Monday 25 January 2010

Linear Discriminant Analysis in C#