判别分析是可用于分类和降维的方法。

由Kaizong Ye，Sherry Deng撰写

线性判别分析（LDA）特别受欢迎，因为它既是分类器又是降维技术。二次判别分析（QDA）是LDA的变体，允许数据的非线性分离。最后，正则化判别分析（RDA）是LDA和QDA之间的折衷。

可下载资源

本文主要关注LDA，并探讨其在理论和实践中作为分类和可视化技术的用途。由于QDA和RDA是相关技术，我不久将描述它们的主要属性以及如何在R中使用它们。

1 文本生成

1.1 文本生成过程

对于一篇文档，可以看成是一组有序的词的序列 $d=(\omega_1,\omega_2,…,\omega_n)$ 。从统计学的角度看，文档的生成可以看成是上帝抛掷骰子生成的结果，每一次抛掷骰子都生成一个单词，抛掷 $n$ 次生成一篇文档。对于该抛掷的过程，存在两个不同学派的观点：

频率学派认为，上帝只有一个骰子，这个骰子有V面，每个面对应一个单词。上帝生成文档是通过独立抛掷 $n$ 次产生n个词；
贝叶斯学派认为，上帝有一个装有多个骰子的坛子，每个骰子有V面，每个面对应一个单词（不同骰子每个面概率不同）。上帝生成文档是通过先在坛子里拿一个骰子，然后再独立抛掷 $n$ 次产生 $n$ 个词。

LDA（Latent Dirichlet Allocation）是Blei等人基于上述贝叶斯学派的观点提出的一种生成文档的方法，该方法假设每个单词是由背后的一个潜在隐藏的主题中抽取出来。对于语料库中的每篇文档，其生成过程为：（1）选择主题；（2）生成单词；（3）组成文档。具体过程如下：

对于每篇文档，其主题存在一个分布，根据该主题分布抽取一个主题；
对于每个主题，其单词存在一个分布，根据该分布抽取一个单词；
重复过程2直到遍历文档中的每一个词。

那么上述主题分布和单词分布都是什么分布呢？且分布需要满足什么性质呢？

1.2 贝叶斯规则

对于一篇文档，构成文档的单词是可观测的，但文档的主题是不可观测的，因此我们需要根据可观测的单词去估计隐藏的主题分布。根据贝叶斯规则：

$p(\theta|x)=\frac{p(x|\theta)*p(\theta)}{p(x)}=\frac{p(x|\theta)p(\theta)}{\int p(x|\theta)*p(\theta)d\theta}\propto p(x|\theta)*p(\theta) \\$

其中：

$p(\theta|x)$ 表示根据观测的数据，得出主题的后验分布；
$p(\theta)$ 表示主题的先验分布；
$p(x|\theta)$ 表示观测数据的似然函数；

上述就是贝叶斯思想：对于一个未知的分布，先假设其服从先验分布（来源于以前做试验数据计算得到，或来自于人们的主观经验），通过观测到的数据，根据贝叶斯规则计算对应的后验分布。很显然 $t$ 时刻的后验分布为 $t+1$ 次的先验分布，因此需要满足先验分布和后验分布具有相同的形式，即上述中 $p(\theta|x)$ 与 $p(\theta)$ 具有相同的分布形式。

2 前置知识

2.1 gamma函数

根据阶乘可知： $n!=1*2*3*…*n$ ，但是当出现小数时，比如计算 $0.5!$ 则阶乘无法直接计算，因此才出现了阶乘的函数形式——gamma函数：

$\Gamma(x)\int_0^{+\infty}e^{-t}t^{x-1}dt(x>0) \\”> </p> <p> 根据gamma函数计算可得<img decoding=$ ，因此有 $\Gamma(n)=(n-1)!$

2.2 二项分布

对于 $n$ 次独立试验，假设每次试验的结果只有两种：成功与失败。成功的概率为 $p$ ，则失败的概率为 $1-p$ 。 $n$ 次独立实验的结果中，成功的次数为 $n_1$ ，失败的次数为 $n_2$ ，则定义似然函数为：

$L=\frac{n!}{n_1!n_2!}p^{n_1}（1-p）^{n_2} \\$

似然函数表示观测数据出现的概率，上式满足： $n=n_1+n_2$

2.3 beta分布

$\begin{aligned} P(p|\alpha_1,\alpha_2)&=\frac{p^{\alpha_1-1}(1-p)^{\alpha_2-1}}{B(\alpha_1,\alpha_2)}\\ &=\frac{\Gamma(\alpha_1+\alpha_2)}{\Gamma(\alpha_1)\Gamma(\alpha_2)}p^{\alpha_1-1}(1-p)^{\alpha_2-1}\\ \end{aligned} \\$

其中 $\alpha_1,\alpha_2$ 为beta分布的超参数，表示伪计数。

参数解释：假设NBA运动员的投篮命中率服从beta分布，对一个刚进联盟的篮球运动员的命中率进行建模，很明显该运动员的命中率也服从beta分布。此时利用上赛季全联盟所有人的投篮命中数 $\alpha_1$ ，打铁数 $\alpha_2$ ，作为该运动员命中率beta分布的参数。且以 $\frac{\alpha_1}{\alpha_1+\alpha_2}$ 作为该运动员的先验命中率。此后，根据该运动员投篮数据的增加，去更新其投篮命中率的计算。

期望计算：

$\begin{aligned} E(p)&=\int_0^1t*Beta(t|\alpha_1,\alpha_2)dt\\ &=\int_0^1t\frac{\Gamma(\alpha_1+\alpha_2)}{\Gamma(\alpha_1)\Gamma(\alpha_2)}t^{\alpha_1-1}(1-t)^{\alpha_2-1}dt\\ &=\frac{\Gamma(\alpha_1+\alpha_2)}{\Gamma(\alpha_1)\Gamma(\alpha_2)}\int_0^1t^{\alpha_1}(1-t)^{\alpha_2-1}dt\\ \end{aligned} \\$

对于分布 $Beta(t|\alpha_1+1,\alpha_2)$ ，有

$\int_0^1Beta(t|\alpha_1+1,\alpha_2)dt=\int_0^1\frac{\Gamma(\alpha_1+1+\alpha_2)}{\Gamma(\alpha_1+1)\Gamma(\alpha_2)}t^{\alpha_1}(1-t)^{\alpha_2-1}dt=1 \\$

则有

$\int_0^1t^{\alpha_1}(1-t)^{\alpha_2-1}dt=\frac{\Gamma(\alpha_1+1)\Gamma(\alpha_2)}{\Gamma(\alpha_1+1+\alpha_2)} \\$

代入上式 $E(p)$ 得

$E(p)=\frac{\Gamma(\alpha_1+\alpha_2)}{\Gamma(\alpha_1)\Gamma(\alpha_2)}*\frac{\Gamma(\alpha_1+1)\Gamma(\alpha_2)}{\Gamma(\alpha_1+1+\alpha_2)}=\frac{\alpha_1}{\alpha_1+\alpha_2} \\$

2.4 共轭先验分布

上述二项式分布的似然函数为

$L=\frac{n!}{n_1!n_2!}p^{n_1}（1-p）^{n_2} \\$

先验的beta分布为

$P(p|\alpha_1,\alpha_2)=\frac{p^{\alpha_1-1}(1-p)^{\alpha_2-1}}{B(\alpha_1,\alpha_2)} \\$

根据贝叶斯规则可得后验分布

$\begin{aligned} P(p|n_1,n_2,\alpha_1,\alpha_2)&=\frac{(\frac{p^{\alpha_1-1}(1-p)^{\alpha_2-1}}{B(\alpha_1,\alpha_2)})*(\frac{n!}{n_1!n_2!}p^{n_1}（1-p）^{n_2})}{\int_0^1(\frac{q^{\alpha_1-1}(1-q)^{\alpha_2-1}}{B(\alpha_1,\alpha_2)})*(\frac{n!}{n_1!n_2!}q^{n_1}（1-q）^{n_2})dq}\\ &=\frac{p^{n_1+\alpha_1-1}(1-p)^{n_2+\alpha_2-1}}{\int_0^1q^{n_1+\alpha_1-1}(1-q)^{n_2+\alpha_2-1}dq}\\ &=\frac{p^{n_1+\alpha_1-1}(1-p)^{n_2+\alpha_2-1}}{B(n_1+\alpha_1,n_2+\alpha_2)}\\ &=Beta(n_1+\alpha_1,n_2+\alpha_2) \end{aligned} \\$

即后验分布仍满足beta分布，与先验一致。

先验分布

$P(p|\alpha_1,\alpha_2)=Beta(\alpha_1,\alpha_2) \\$

后验分布

$P(p|\alpha_1,\alpha_2,n_1,n_2)=Beta(\alpha_1+n_1,\alpha_2+n_2) \\$

其中 $\alpha_1,\alpha_2$ 为伪计数， $n_1,n_2$ 为观测的数据。

综上，二项分布对应的共轭先验为beta分布，其意思为，beta分布*二项分布的结果仍未beta分布。

3 Dirichlet分布

3.1 多项式分布

对二项式分布推广到k种结果的情况，此时变为多项式分布，对应的似然函数为

$L=\frac{n!}{n_1!n_2!…n_k!}p_1^{n_1}p_2^{n_2}…p_k^{n_k}=n!\prod_{i=1}^k\frac{p_i^{n_i}}{n_i!} \\$

其中， $\sum_{i=1}^kn_i=n$ ， $\sum_{i=1}^kp_i=1$

3.2 Dirichlet分布

对beta分布推广到k种结果的情况，此时表为Dirichle分布

$P(p_1,p_2,…,p_k|\alpha_1,\alpha_2,…,\alpha_k)=\frac{\Gamma(\alpha_1+\alpha_2+…+\alpha_k)}{\Gamma(\alpha_1)\Gamma(\alpha_2)…\Gamma(\alpha_k)}p_1^{\alpha_1-1}p_2^{\alpha_2-1}…p_k^{\alpha_k-1} \\$

其中 $\sum_{i=1}^kp_i=1$

同2.3，Dirichlet分布的期望为： $E(p_i)=\frac{\alpha_i}{\sum_{i=1}^k\alpha_i}$ 。

3.3 多项式分布与Dirichlet分布共轭

同2.4可得

先验分布

$P(p_1,p_2,…,p_k|\alpha_1,\alpha_2,…,\alpha_k)=Beta(\alpha_1,\alpha_2,…,\alpha_k) \\$

后验分布

$P(p_1,p_2,…,p_k|\alpha_1,\alpha_2,…,\alpha_k,n_1,n_2,…,n_k)=Beta(\alpha_1+n_1,\alpha_2+n_2,…,\alpha_k+n_k) \\$

其中 $\alpha_1,\alpha_2,…,\alpha_k$ 为伪计数， $n_1,n_2,…,n_k$ 为观测的数据。

4 LDA模型

4.1 模型表示

现在问题是这样的，我们有 $D$ 篇文档，对应第 $d$ 个文档有 $N_d$ 个单词。我们的目标是找到每一篇文档的主题分布和每一个主题单词的分布。在LDA模型中，我们需要先假定一个主题数目 $K$ ，这样所有的分布就都基于 $K$ 个主题展开。

对于语料库中的每篇文档，LDA定义了一个生成过程，以1文本生成部分的投骰子为例，如下：

$\vec\alpha\to\vec\theta_d\to z_{d,n}$ 。这个过程表示在生成第 $d$ 篇文档的时候，在D个服从 $\vec\theta_d～Dirichlet(\vec\alpha)$ 的坛子中找到骰子 $\vec\theta_d$ ，然后投掷这枚骰子生成第n个词的主题（topic）编号 $z_{d,n}$ ；
$\vec\eta\to\beta_k\to\omega_{d,n}|k=z_{d,n}$ 。这个过程表示在生成第 $d$ 篇文档第 $n$ 个词的时候，在K个服从 $\vec\beta_k～Dirichlet(\vec\eta)$ 的坛子找到主题编号为 $z_{d,n}$ 的骰子，然后投掷这枚骰子生成单词 $w_{d,n}$

简单来说，步骤1就是 $p(\vec z)$ ，步骤2就是 $p(\vec w|\vec z)$ ，根据条件概率的基本公式可得：

$p(\vec w,\vec z)=p(\vec w|\vec z)p(\vec z) \\$

而LDA的目标就是要找出每个词后潜在的主题，所以为了达到这个目标，需要计算后验概率：

$p(\vec z|\vec w)=\frac{p(\vec w,\vec z)}{\sum_zp(\vec w,\vec z)} \\$

4.2 直接计算

针对上述后验概率的计算，在这里对其计算复杂度进行分析。

按照离散分布边缘概率的处理方式，文档中一个单词 $w_i$ 的全概率为

$p(w_i)=\sum_{k=1}^Kp(w_i,z_i=k)=\sum_{k=1}^Kp(w_i|z_i=k)p(z_i=k) \\$

因此对于上述后验概率计算中的分母 $\sum_zp(\vec w,\vec z)$ ，其表示为所有单词的联合概率，则有

$\sum_zp(\vec w,\vec z)=p(\vec w)=\prod_{i=1}^Vp(w_i)=\prod_{i=1}^V\sum_{k=1}^Kp(w_i|z_i=k)p(z_i=k) \\$

每个单词都对应着 $K$ 个主题，总共有 $V$ 个单词，因此分母计算陷入了 $K^V$ 项难题，这个离散状态空间太大以至于无法枚举。

5 模型求解

在上节结尾可知，现实中，往往很难求出精确的概率，因此常常采用近似推断的方法。近似推断的方法大致可分为两大类：第一类是使用随机化采样完成近似，比如Gibbs采样；第二类是使用确定性近似完成近似推断，比如变分推断。

线性判别分析

LDA是一种分类和降维技术，可以从两个角度进行解释。第一个是解释是概率性的，第二个是更多的程序解释，归功于费舍尔。第一种解释对于理解LDA的假设是有用的。第二种解释可以更好地理解LDA如何降低维数。

Fisher的LDA优化标准

Fisher的LDA优化标准规定组的质心应尽可能分散。这相当于找到一个线性组合ž= aŤXZ=aTX，使得aTaT相对于类内方差的类间方差最大化。

LDA模型的复杂性

LDA的有效参数的数量可以通过以下方式导出。有KK手段μ^kμ^k被估计。协方差矩阵不需要额外的参数，因为它已经由质心定义。由于我们需要估计KK判别函数（以获得判定边界），这就产生了涉及p个元素的KK计算。另外，我们有ķ-1为自由参数ķ前科。因此，有效LDA参数的数量是Kp+（K-1）。

课程

R语言数据分析挖掘必知必会

从数据获取和清理开始，有目的的进行探索性分析与可视化。让数据从生涩的资料，摇身成为有温度的故事。

立即参加 ➜

LDA摘要

在这里，我总结了LDA的两个观点，并总结了该模型的主要特性。

概率论

LDA使用贝叶斯规则来确定观察xx属于kk类的后验概率。由于LDA的正常假设，后验由多元高斯定义，其协方差矩阵假定对于所有类是相同的。新的点通过计算判别函数分类δkδk（后验概率的枚举器）并返回类kk具有最大δkδk。判别变量可以通过类内和类间方差的特征分解来获得。

费舍尔的观点

根据Fisher，LDA可以理解为降维技术，其中每个连续变换是正交的并且相对于类内方差最大化类间方差。此过程将特征空间转换为具有K−1维度的仿射空间。在对输入数据进行扩展之后，可以通过在考虑类先验的情况下确定仿射空间中的最接近的质心来对新点进行分类。

LDA的特性

LDA具有以下属性：

LDA假设数据是高斯数据。更具体地说，它假定所有类共享相同的协方差矩阵。
LDA在K−1维子空间中找到线性决策边界。因此，如果自变量之间存在高阶相互作用，则不适合。
LDA非常适合于多类问题，但是当类分布不平衡时应该小心使用，因为根据观察到的计数来估计先验。因此，观察很少被分类为不常见的类别。
与PCA类似，LDA可用作降维技术。请注意，LDA的转换本质上与PCA不同，因为LDA是一种考虑结果的监督方法。

数据集

为了举例说明线性判别分析，我们将使用音素语音识别数据集。该数据集对于展示判别分析很有用，因为它涉及五种不同的结果。

library(RCurl)
f <- getURL('phoneme.csv')
df <- read.csv(textConnection(f), header=T)
print(dim(df))

## [1] 4509  259

为了以后评估模型，我们将每个样本分配到训练或测试集中：

最受欢迎的见解

1.用SPSS估计HLM层次线性模型模型

2.R语言线性判别分析（LDA），二次判别分析（QDA）和正则判别分析（RDA）

3.基于R语言的lmer混合线性回归模型

4.R语言Gibbs抽样的贝叶斯简单线性回归仿真分析

5.在r语言中使用GAM（广义相加模型）进行电力负荷时间序列分析

6.使用SAS，Stata，HLM，R，SPSS和Mplus的分层线性模型HLM

7.R语言中的岭回归、套索回归、主成分回归：线性模型选择和正则化

8.R语言用线性回归模型预测空气质量臭氧数据

9.R语言分层线性模型案例

#逻辑变量，判断是否是训练集
train <- grepl("^train", df$speaker)
# 移除非特征变量
to.exclude <- c("row.names", "speaker""g")
feature.df <- df[, !colnames(df) %<strong>in</strong>% to.exclude]
test.set <- subset(feature.df, !train)
train.set <- subset(feature.df, train)
train.responses <- subset(df, train)$g
test.responses <- subset(df, !train)$g

在R中拟合LDA模型

我们可以通过以下方式拟合LDA模型：

 lda.model <- lda(train.set, grouping = train.responses)

R语言主成分回归（PCR）、多元线性回归特征降维分析光谱数据和汽车油耗、性能数据

阅读文章 ➜

LDA作为可视化技术

我们可以通过在缩放数据上应用变换矩阵将训练数据转换为规范坐标。要获得与predict.lda函数返回的结果相同的结果，我们需要首先围绕加权平均数据居中：

## [1] TRUE

我们可以使用前两个判别变量来可视化数据：

绘制两个LDA维度中的数据显示三个集群：

群集1（左）由aa和ao音素组成
群集2（右下角）由dcl和iy音素组成
群集3（右上角）由sh音素组成

这表明两个维度不足以区分所有5个类别。然而，聚类表明可以非常好地区分彼此充分不同的音素。

我们还可以使用plot.lda函数绘制训练数据到所有判别变量对的映射，其中dimen参数可用于指定所考虑的维数：

随时关注您喜欢的主题

为了可视化组的质心，我们可以创建自定义图：

解释后验概率

除了将数据转换为由分量x提供的判别变量之外，预测函数还给出后验概率，其可以用于分类器的进一步解释。例如：

## [1] "Posterior of predicted class 'sh' is: 1"
##        aa    ao   dcl    iy    sh
## aa  0.797 0.203 0.000 0.000 0.000
## ao  0.123 0.877 0.000 0.000 0.000
## dcl 0.000 0.000 0.985 0.014 0.002
## iy  0.000 0.000 0.001 0.999 0.000
## sh  0.000 0.000 0.000 0.000 1.000

各个班级的后验表格表明该模型对音素aa和ao最不确定，这与我们对可视化的期望一致。

LDA作为分类器

如前所述，LDA的好处是我们可以选择用于分类的规范变量的数量。在这里，我们仍将通过使用多达四个规范变量进行分类来展示降级LDA的使用。

##   Rank Accuracy
## 1    1     0.51
## 2    2     0.71
## 3    3     0.86
## 4    4     0.92

由于具有四维的LDA获得最大精度，我们将决定使用所有判别坐标进行分类。

正如从变换空间的视觉探索所预期的那样，测试精度随着每个附加维度而增加。

为了解释模型，我们可以可视化分类器的性能：

在图中，预测的音素以不同的颜色显示，而模型预测通过不同的符号显示。具有100％准确度的模型将为每种颜色分配单个符号。

二次判别分析

QDA是LDA的变体，其中针对每类观察估计单个协方差矩阵。如果事先知道个别类别表现出不同的协方差，则QDA特别有用。QDA的缺点是它不能用作降维技术。

由于QDA估计每个类的协方差矩阵，因此它具有比LDA更多的有效参数。我们可以通过以下方式得出参数的数量。

因此，QDA参数的有效数量是ķ- 1 + K.p + K.p （p + 1 ）2K−1+Kp+Kp(p+1)2。

由于QDA参数的数量在pp是二次的，因此当特征空间很大时，应小心使用QDA。

QDA在R

我们可以通过以下方式执行QDA：

的QDA和LDA对象之间的主要区别是，QDA具有p×p的变换矩阵对于每个类k∈{1,…,K}。这些矩阵确保组内协方差矩阵是球形的，但不会导致子空间减小。因此，QDA不能用作可视化技术。

## [1] "Accuracy of QDA is: 0.84"

QDA的准确度略低于全级LDA的准确度。这可能表明共同协方差的假设适合于该数据集。

让我们确定QDA在音素数据集上是否优于LDA：

规范的判别分析

由于RDA是一种正则化技术，因此当存在许多潜在相关的特征时。现在让我们评估音素数据集上的RDA。

R中的RDA

rda.preds <- predict(rda.model, t(train.set), train.responses, t(test.set))

# 判断alpha的表现
rda.perf <- vector(, dim(rda.preds)[1])
for(i in seq(dim(rda.preds)[1])) {
     # 判断gamma的表现
    res <- apply(rda.preds[i,,], 1, function(x) length(which(x == as.numeric(test.responses))) / length(test.responses))
    rda.perf[[i]] <- res
}
rda.perf <- do.call(rbind, rda.perf)
rownames(rda.perf) <- alphas

结论

判别分析对于多类问题特别有用。LDA非常易于理解，因为它可以减少维数。使用QDA，可以建模非线性关系。RDA是一种正则化判别分析技术，对大量特征特别有用。

可下载资源

关于作者

Kaizong Ye是拓端研究室（TRL）的研究员。在此对他对本文所作的贡献表示诚挚感谢，他在上海财经大学完成了统计学专业的硕士学位，专注人工智能领域。擅长Python.Matlab仿真、视觉处理、神经网络、数据分析。

本文借鉴了作者最近为《R语言数据分析挖掘必知必会 》课堂做的准备。

非常感谢您阅读本文，如需帮助请联系我们！

R语言线性判别分析（LDA），二次判别分析（QDA）和正则判别分析（RDA）