What is the relationship between regression and linear discriminant analysis (LDA)? Is there a relationship between regression and linear discriminant analysis (LDA)? What are their similarities and differences? Does it make any difference if there are two classes or more than two classes?
 A: The purpose of this answer is to explain the exact mathematical relationship between linear discriminant analysis (LDA) and multivariate linear regression (MLR). It will turn out that the correct framework is provided by reduced rank regression (RRR).
We will show that LDA is equivalent to RRR of the whitened class indicator matrix on the data matrix.

Notation
Let $\newcommand{\X}{\mathbf X}\X$ be the $n\times d$ matrix with data points $\newcommand{\x}{\mathbf x}\x_i$ in rows and variables in columns. Each point belongs to one of the $k$ classes, or groups. Point $\x_i$ belongs to class number $g(i)$.
Let $\newcommand{\G}{\mathbf G}\G$ be the $n \times k$ indicator matrix encoding group membership as follows: $G_{ij}=1$ if $\x_i$ belongs to class $j$, and $G_{ij}=0$ otherwise. There are $n_j$ data points in class $j$; of course $\sum n_j = n$.
We assume that the data are centered and so the global mean is equal to zero, $\newcommand{\bmu}{\boldsymbol \mu}\bmu=0$. Let $\bmu_j$ be the mean of class $j$.
LDA
The total scatter matrix $\newcommand{\C}{\mathbf C}\C=\X^\top \X$ can be decomposed into the sum of between-class and within-class scatter matrices defined as follows:
\begin{align}
\C_b &= \sum_j n_j \bmu_j \bmu_j^\top \\
\C_w &= \sum(\x_i - \bmu_{g(i)})(\x_i - \bmu_{g(i)})^\top.
\end{align}
One can verify that $\C = \C_b + \C_w$. LDA searches for discriminant axes that have maximal between-group variance and minimal within-group variance of the projection. Specifically, first discriminant axis is the unit vector $\newcommand{\w}{\mathbf w}\w$ maximizing $\w^\top \C_b \w / (\w^\top \C_w \w)$, and the first $p$ discriminant axes stacked together into a matrix $\newcommand{\W}{\mathbf W}\W$ should maximize the trace $$\DeclareMathOperator{\tr}{tr} L_\mathrm{LDA}=\tr\left(\W^\top \C_b \W (\W^\top \C_w \W)^{-1}\right).$$
Assuming that $\C_w$ is full rank, LDA solution $\W_\mathrm{LDA}$ is the matrix of eigenvectors of $\C_w^{-1} \C_b$ (ordered by the eigenvalues in the decreasing order).
This was the usual story. Now let us make two important observations.
First, within-class scatter matrix can be replaced by the total scatter matrix (ultimately because maximizing $b/w$ is equivalent to maximizing $b/(b+w)$), and indeed, it is easy to see that $\C^{-1} \C_b$ has the same eigenvectors.
Second, the between-class scatter matrix can be expressed via the group membership matrix defined above. Indeed, $\G^\top \X$ is the matrix of group sums. To get the matrix of group means, it should be multiplied by a diagonal matrix with $n_j$ on the diagonal; it's given by $\G^\top \G$. Hence, the matrix of group means is $(\G^\top \G)^{-1}\G^\top \X$ (sapienti will notice that it's a regression formula). To get $\C_b$ we need to take its scatter matrix, weighted by the same diagonal matrix, obtaining $$\C_b = \X^\top \G (\G^\top \G)^{-1}\G^\top \X.$$  If all $n_j$ are identical and equal to $m$ ("balanced dataset"), then this expression simplifies to $\X^\top \G \G^\top \X / m$.
We can define normalized indicator matrix $\newcommand{\tG}{\widetilde {\mathbf G}}\tG$ as having $1/\sqrt{n_j}$ where $\G$ has $1$. Then for both, balanced and unbalanced datasets, the expression is simply $\C_b = \X^\top \tG \tG^\top \X$. Note that $\tG$ is, up to a constant factor, the whitened indicator matrix: $\tG = \G(\G^\top \G)^{-1/2}$.
Regression
For simplicity, we will start with the case of a balanced dataset. 
Consider linear regression of $\G$ on $\X$. It finds $\newcommand{\B}{\mathbf B}\B$ minimizing $\| \G - \X \B\|^2$. Reduced rank regression does the same under the constraint that $\B$ should be of the given rank $p$. If so, then $\B$ can be written as $\newcommand{\D}{\mathbf D} \newcommand{\F}{\mathbf F} \B=\D\F^\top$ with both $\D$ and $\F$ having $p$ columns. One can show that the rank two solution can be obtained from the rank solution by keeping the first column and adding an extra column, etc.
To establish the connection between LDA and linear regression, we will prove that $\D$ coincides with $\W_\mathrm{LDA}$.
The proof is straightforward. For the given $\D$, optimal $\F$ can be found via regression: $\F^\top = (\D^\top \X^\top \X \D)^{-1} \D^\top \X^\top \G$. Plugging this into the loss function, we get $$\| \G - \X \D (\D^\top \X^\top \X \D)^{-1} \D^\top \X^\top \G\|^2,$$ which can be written as trace using the identity $\|\mathbf A\|^2=\mathrm{tr}(\mathbf A \mathbf A^\top)$. After easy manipulations we get that the regression is equivalent to maximizing (!) the following scary trace: $$\tr\left(\D^\top \X^\top \G \G^\top \X \D (\D^\top \X^\top \X \D)^{-1}\right),$$ which is actually nothing else than $$\ldots = \tr\left(\D^\top \C_b \D (\D^\top \C \D)^{-1}\right)/m \sim L_\mathrm{LDA}.$$
This finishes the proof. For unbalanced datasets we need to replace $\G$ with $\tG$.
One can similarly show that adding ridge regularization to the reduced rank regression is equivalent to the regularized LDA.
Relationship between LDA, CCA, and RRR
In his answer, @ttnphns made a connection to canonical correlation analysis (CCA). Indeed, LDA can be shown to be equivalent to CCA between $\X$ and $\G$. In addition, CCA between any $\newcommand{\Y}{\mathbf Y}\Y$ and $\X$ can be written as RRR predicting whitened $\Y$ from $\X$. The rest follows from this. 
Bibliography
It is hard to say who deserves the credit for what is presented above.
There is a recent conference paper by Cai et al. (2013) On The Equivalent of Low-Rank Regressions and Linear Discriminant Analysis Based Regressions that presents exactly the same proof as above but creates the impression that they invented this approach. This is definitely not the case. Torre wrote a detailed treatment of how most of the common linear multivariate methods can be seen as reduced rank regression, see A Least-Squares Framework for Component Analysis, 2009, and a later book chapter A unification of component analysis methods, 2013; he presents the same argument but does not give any references either. This material is also covered in the textbook Modern Multivariate Statistical Techniques (2008) by Izenman, who introduced RRR back in 1975.
The relationship between LDA and CCA apparently goes back to Bartlett, 1938, Further aspects of the theory of multiple regression -- that's the reference I often encounter (but did not verify). The relationship between CCA and RRR is described in the Izenman, 1975, Reduced-rank regression for the multivariate linear model. So all of these ideas have been around for a while.
A: I take it that the question is about LDA and linear (not logistic) regression.
There is a considerable and meaningful relation between linear regression and linear discriminant analysis. In case the dependent variable (DV) consists just of 2 groups the two analyses are actually identical. Despite that computations are different and the results - regression and discriminant coefficients - are not the same, they are exactly proportional to each other.
Now for the more-than-two-groups situation. First, let us state that LDA (its extraction, not classification stage) is equivalent (linearly related results) to canonical correlation analysis if you turn the grouping DV into a set of dummy variables (with one redundant of them dropped out) and do canonical analysis with sets "IVs" and "dummies". Canonical variates on the side of "IVs" set that you obtain are what LDA calls "discriminant functions" or "discriminants".
So, then how canonical analysis is related to linear regression? Canonical analysis is in essence a MANOVA (in the sense "Multivariate Multiple linear regression" or "Multivariate general linear model") deepened into latent structure of relationships between the DVs and the IVs. These two variations are decomposed in their inter-relations into latent "canonical variates". Let us take the simplest example, Y vs X1 X2 X3. Maximization of correlation between the two sides is linear regression (if you predict Y by Xs) or - which is the same thing - is MANOVA (if you predict Xs by Y). The correlation is unidimensional (with magnitude R^2 = Pillai's trace) because the lesser set, Y, consists just of one variable. Now let's take these two sets: Y1 Y2 vs X1 x2 x3. The correlation being maximized here is 2-dimensional because the lesser set contains 2 variables. The first and stronger latent dimension of the correlation is called the 1st canonical correlation, and the remaining part, orthogonal to it, the 2nd canonical correlation. So, MANOVA (or linear regression) just asks what are partial roles (the coefficients) of variables in the whole 2-dimensional correlation of sets; while canonical analysis just goes below to ask what are partial roles of variables in the 1st correlational dimension, and in the 2nd.
Thus, canonical correlation analysis is multivariate linear regression deepened into latent structure of relationship between the DVs and IVs. Discriminant analysis is a particular case of canonical correlation analysis (see exactly how). So, here was the answer about the relation of LDA to linear regression in a general case of more-than-two-groups.
Note that my answer does not at all see LDA as classification technique. I was discussing LDA only as extraction-of-latents technique. Classification is the second and stand-alone stage of LDA (I described it here). @Michael Chernick was focusing on it in his answers.
A: Here is a reference to one of Efron's papers: The Efficiency of Logistic Regression Compared to Normal Discriminant Analysis, 1975.
Another relevant paper is Ng & Jordan, 2001, On Discriminative vs. Generative classifierers: A comparison of logistic regression and naive Bayes. And here is an abstract of a comment on it by Xue & Titterington, 2008, that mentions O'Neill's papers related to his PhD dissertation:

Comparison of generative and discriminative classifiers is an
  ever-lasting topic. As an important contribution to this topic, based
  on their theoretical and empirical comparisons between the naïve Bayes
  classifier and linear logistic regression, Ng and Jordan (NIPS
  841---848, 2001) claimed that there exist two distinct regimes of
  performance between the generative and discriminative classifiers with
  regard to the training-set size. In this paper, our empirical and
  simulation studies, as a complement of their work, however, suggest
  that the existence of the two distinct regimes may not be so reliable.
  In addition, for real world datasets, so far there is no theoretically
  correct, general criterion for choosing between the discriminative and
  the generative approaches to classification of an observation $x$ into
  a class $y$; the choice depends on the relative confidence we have in
  the correctness of the specification of either $p(y|x)$ or $p(x, y)$
  for the data. This can be to some extent a demonstration of why Efron
  (J Am Stat Assoc 70(352):892---898, 1975) and O'Neill (J Am Stat Assoc
  75(369):154---160, 1980) prefer normal-based linear discriminant
  analysis (LDA) when no model mis-specification occurs but other
  empirical studies may prefer linear logistic regression instead.
  Furthermore, we suggest that pairing of either LDA assuming a common
  diagonal covariance matrix (LDA) or the naïve Bayes classifier and
  linear logistic regression may not be perfect, and hence it may not be
  reliable for any claim that was derived from the comparison between
  LDA or the naïve Bayes classifier and linear logistic regression to be
  generalised to all generative and discriminative classifiers.

There are a lot of other references on this that you can find online. 
A: Linear regression and linear discriminant analysis are very different.  Linear regression relates a dependent variable to a set of independent predictor variables.  The idea is to find a function linear in the parameters that best fits the data.  It does not even have to be linear in the covariates.  Linear discriminant analysis on the other hand is a procedure for classifying objects into categories. For the two-class problem it seeks to find the best separating hyperplane for dividing the groups into two catgories.  Here best means that it minimizes a loss function that is a linear combination of the error rates.  For three or more groups it finds the best set of hyperplanes (k-1 for the k class problem).  In discriminant analysis the hypoerplanes are linear in the feature variables.
The main similarity between the two is term linear in the titles.
