19
$\begingroup$

In R, I use lda function from library MASS to do classification. As I understand LDA, input $x$ will be assigned label $y$, which maximize $p(y|x)$, right?

But when I fit the model, in which $$x=(Lag1,Lag2)$$$$y=Direction,$$ I don't quite understand the output from lda,

Edit: to reproduce the output below, first run:

library(MASS)
library(ISLR)

train = subset(Smarket, Year < 2005)

lda.fit = lda(Direction ~ Lag1 + Lag2, data = train)
> lda.fit
Call:
lda(Direction ~ Lag1 + Lag2, data = train)

Prior probabilities of groups:
    Down       Up 
0.491984 0.508016 

Group means:
            Lag1        Lag2
Down  0.04279022  0.03389409
Up   -0.03954635 -0.03132544

Coefficients of linear discriminants:
            LD1
Lag1 -0.6420190
Lag2 -0.5135293

I understand all the info in the above output but one thing, what is LD1? I search the web for it, is it linear discriminant score? What is that and why do I need it?

UPDATE

I read several posts (such as this and this one) and also search the web for DA, and now here is what I think about DA or LDA.

  1. It can be used to do classification, and when this is the purpose, I can use the Bayes approach, that is, compute the posterior $p(y|x)$ for each class $y_i$, and then classify $x$ to the class with the highest posterior. By this approach, I don't need to find out the discriminants at all, right?

  2. As I read in the posts, DA or at least LDA is primarily aimed at dimensionality reduction, for $K$ classes and $D$-dim predictor space, I can project the $D$-dim $x$ into a new $(K-1)$-dim feature space $z$, that is, \begin{align*}x&=(x_1,...,x_D)\\z&=(z_1,...,z_{K-1})\\z_i&=w_i^Tx\end{align*}, $z$ can be seen as the transformed feature vector from the original $x$, and each $w_i$ is the vector on which $x$ is projected.

Am I right about the above statements? If yes, I have following questions:

  1. What is a discriminant? Is each entry $z_i$ in vector $z$ is a discriminant? Or $w_i$?

  2. How to do classification using discriminants?

$\endgroup$
13
  • 1
    $\begingroup$ LDA has 2 distinct stages: extraction and classification. At extraction, latent variables called discriminants are formed, as linear combinations of the input variables. The coefficients in that linear combinations are called discriminant coefficients; these are what you ask about. On the 2nd stage, data points are assigned to classes by those discriminants, not by original variables. To read more, search discriminant analysis on this site. $\endgroup$
    – ttnphns
    Commented Feb 22, 2014 at 7:51
  • 2
    $\begingroup$ Linear discriminant score is a value of a data point by a discriminant, so don't confuse it with discriminant coefficient, which is like a regressional coefficient. See my detailed answer here. $\endgroup$
    – ttnphns
    Commented Feb 22, 2014 at 8:09
  • $\begingroup$ @ttnphns, thanks and I'll read more about DA. BTW, I thought that to classify an input $X$, I just need to compute the posterior $p(y|x)$ for all the classes and then pick the class with highest posterior, right? And I don't see why I need $LD1$ in the computation of posterior. $\endgroup$
    – avocado
    Commented Feb 22, 2014 at 8:47
  • $\begingroup$ You can and may do Bayes-rule classification based on the original variables. But this won't be discriminant analysis. The essential part of LDA is that dimensionality reduction, which allows you to replace the original variables-classifiers by a smaller number of derivative classifiers, the discriminants. Please read posts here, particularly mine, they thorougly describe ideas and maths of LDA. $\endgroup$
    – ttnphns
    Commented Feb 22, 2014 at 9:16
  • $\begingroup$ @ttnphns, I'm reading the post you linked in the above comment, ;-) $\endgroup$
    – avocado
    Commented Feb 22, 2014 at 10:14

3 Answers 3

7
$\begingroup$

If you multiply each value of LDA1 (the first linear discriminant) by the corresponding elements of the predictor variables and sum them ($-0.6420190\times$Lag1$+ -0.5135293\times$Lag2) you get a score for each respondent. This score along the the prior are used to compute the posterior probability of class membership (there are a number of different formulas for this). Classification is made based on the posterior probability, with observations predicted to be in the class for which they have the highest probability.

The chart below illustrates the relationship between the score, the posterior probability, and the classification, for the data set used in the question. The basic patterns always holds with two-group LDA: there is 1-to-1 mapping between the scores and the posterior probability, and predictions are equivalent when made from either the posterior probabilities or the scores.

Score, Posterior Probability, Classification

Answers to the sub-questions and some other comments

  • Although LDA can be used for dimension reduction, this is not what is going on in the example. With two groups, the reason only a single score is required per observation is that this is all that is needed. This is because the probability of being in one group is the complement of the probability of being in the other (i.e., they add to 1). You can see this in the chart: scores of less than -.4 are classified as being in the Down group and higher scores are predicted to be Up.

  • Sometimes the vector of scores is called a discriminant function. Sometimes the coefficients are called this. I'm not clear on whether either is correct. I believe that MASS discriminant refers to the coefficients.

  • The MASS package's lda function produces coefficients in a different way to most other LDA software. The alternative approach computes one set of coefficients for each group and each set of coefficients has an intercept. With the discriminant function (scores) computed using these coefficients, classification is based on the highest score and there is no need to compute posterior probabilities in order to predict the classification. I have put some LDA code in GitHub which is a modification of the MASS function but produces these more convenient coefficients (the package is called Displayr/flipMultivariates, and if you create an object using LDA you can extract the coefficients using obj$original$discriminant.functions).

  • I have posted the R for code all the concepts in this post here.

  • There is no single formula for computing posterior probabilities from the score. The easiest way to understand the options is (for me anyway) to look at the source code, using:

library(MASS) getAnywhere("predict.lda")

$\endgroup$
3
  • $\begingroup$ I'm not clear on whether either [word use] is correct "discriminant function" aka "discriminant" is an extracted variate - a variable, a dimension. It therefore is characterized both by coefficients (weights) to assess it from the input variables, and by scores, the values. Exactly like a PC in PCA. So, "discriminant coefficients" and "discriminant scores" are the correct usage. $\endgroup$
    – ttnphns
    Commented Feb 20, 2018 at 12:16
  • $\begingroup$ @ttnphns, your usage of the terminology is very clear and unambiguous. But, it is not the usage that appears in much of the post and publications on the topic, which is the point that I was trying to make. Based on word-meaning alone, it is pretty clear to me that the "discriminant function" should refer to the mathematical function (i.e., sumproduct and the coefficients), but again it is not clear to me that this is the widespread usage. $\endgroup$
    – Tim
    Commented Feb 20, 2018 at 19:43
  • $\begingroup$ @Tim the link you've posted for the code is dead , can you copy the code into your answer please? $\endgroup$
    – baxx
    Commented Aug 27, 2019 at 12:09
3
$\begingroup$

Discriminant in the context of ISLR, 4.6.3 Linear Discriminant Analysis, pp161-162 is, as I understand, the value of \begin{equation} \hat\delta_2(\vec x) - \hat\delta_1(\vec x) = {\vec x}^T\hat\Sigma^{-1}\Bigl(\vec{\hat\mu}_2 - \vec{\hat\mu}_1\Bigr) - \frac{1}{2}\Bigl(\vec{\hat\mu}_2 + \vec{\hat\mu}_1\Bigr)^T\hat\Sigma^{-1}\Bigl(\vec{\hat\mu}_2 - \vec{\hat\mu}_1\Bigr) + \log\Bigl(\frac{\pi_2}{\pi_1}\Bigr), \tag{$*$} \end{equation}

where $\vec x = (\mathrm{Lag1}, \mathrm{Lag2})^T$. For the 2nd term in $(*)$, it should be noted that, for symmetric matrix M, we have $\vec x^T M\vec y = \vec y^T M \vec x$.

LD1 is the coefficient vector of $\vec x$ from above equation, which is \begin{equation} \hat\Sigma^{-1}\Bigl(\vec{\hat\mu}_2 - \vec{\hat\mu}_1\Bigr). \end{equation} $y$ at $\vec x$ is 2 if $(*)$ is positive, and 1 if $(*)$ is negative.

LD1 is given as lda.fit$scaling. The discriminant vector ${\vec x}^T\hat\Sigma^{-1}\Bigl(\vec{\hat\mu}_2 - \vec{\hat\mu}_1\Bigr)$ computed using LD1 for a test set is given as lda.pred$x, where

test = subset(Smarket, Year==2005)
lda.pred = predict(lda.fit, test)

test set is not necessarily given as above, it can be given arbitrarily.

Unfortunately, lda.pred$x alone cannot tell whether $y$ is 1 or 2. We need the 2nd and the 3rd term in $(*)$. I could not find these terms from the output of lda() and/or predict(lda.fit,..).

We can compute all three terms of $(*)$ by hand, I mean using just the basic functions of R. The script for LD1 is given below.

bTrain = (Smarket$Year < 2005) # boolean vector for the training set
bUp = (Smarket$Direction == "Up") # boolean vector for "Up" group
up.group = Smarket[bTrain & bUp,]
down.group = Smarket[bTrain & !bUp,]
cov.up = cov(up.group[,c(2,3)])
cov.down = cov(down.group[,c(2,3)])
nUp = nrow(up.group)
nDown = nrow(down.group)
n = nUp + nDown; K = 2
Sigma <- 1/(n - K) * (cov.up * (nUp - 1) + cov.down * (nDown - 1))
SigmaInv <- solve(Sigma)
lda.fit = lda(Direction ~ Lag1+Lag2, data=Smarket, subset=bTrain) 
mu2_mu1 = matrix(lda.fit$means[2,] - lda.fit$means[1,])
coeff.vec = SigmaInv %*% mu2_mu1; coeff.vec
            [,1]
Lag1 -0.05544078
Lag2 -0.04434520
v.scalar = sqrt(t(coeff.vec) %*% Sigma %*% coeff.vec)
myLD1 = coeff.vec/drop(v.scalar); myLD1 # same as below
           [,1]
Lag1 -0.6420190
Lag2 -0.5135293
lda.fit$scaling # same as above
            LD1
Lag1 -0.6420190
Lag2 -0.5135293

Here is the catch: myLD1 is perfectly good in the sense that it can be used in classifying $\vec x$ according to the value of its corresponding response variable $y$. Since the discriminant function $(*)$ is linear in $\vec x$ (actually it's not linear, it's affine) any scalar multiple of myLD1 will do the job provided that the second and the third term are multiplied by the same scalar, which is 1/v.scalar in the code above.

LD1 given by lda() has the nice property that the generalized norm is 1, which our myLD1 lacks.

t(lda.fit$scaling) %*% Sigma %*% lda.fit$scaling
    LD1
LD1   1

-end-

$\endgroup$
1
  • $\begingroup$ From formula $(*)$, one can see that the midpoint (mu1 + mu2)/2 lies on the decision boundary in case $\pi_1 = \pi_2$. $\endgroup$ Commented Mar 6, 2020 at 5:17
0
$\begingroup$

The theory behind this function is "Fisher's Method for Discriminating among Several Population". I recommend chapter 11.6 in applied multivariate statistical analysis(ISBN: 9780134995397) for reference.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.