R语言进阶之判别分析

R语言的“MASS“包是一个十分强大的统计包，可以进行各种统计分析，我也将围绕它来介绍判别分析。”MASS“包既可以进行线性判别，也可以进行二次判别。除非指定先验概率，否则”MASS”的判别分析会按照样本量来等比例给出先验概率。

这里我们使用的还是鸢尾花数据集，将“setosa”编号为0，“versicolor”编号为1，“virginica”编号为2。

1. 线性判别函数

# 使用Jacknifed预测进行线性判别分析
library(MASS) # 加载R包
mydata <- iris # 将iris命名成mydata以便后续操作
mydata$type[which(mydata$Species=="setosa")] <- 0 # 设置setosa为0
mydata$type[which(mydata$Species=="versicolor")] <- 1 # 设置versicolor为1
mydata$type[which(mydata$Species=="virginica")] <- 2 # 设置virginica为2
attach(mydata) # 固定数据框
fit <- lda(type~ Sepal.Length+Sepal.Width+Petal.Length+Petal.Width, data=mydata,na.action="na.omit") # 拟合判别模型，忽视缺失值
fit # 展示拟合结果

注意一下，上面的函数lda()其实就是线性判别分析函数（linear discriminant analysis），当指定na.action="na.omit"就相当于删除含有缺失值的样本。

注意：lda()函数并不提供显著性检验，如果想进行显著性检验，可以使用manova()函数。

2. 二次判别函数

在“MASS”包里我们可以使用qda()函数来进行二次判别分析，这里qda实际上就是二次判别分析（quadratic discriminant analysis）的英文缩写。二次判别分析的前提条件比较宽松，不像线性判别那样要求同方差。同样地，我们仍然使用之间建立好的mydata数据框来进行二次判别分析。

# 三组类别的二次判别分析 
# 假定各组的先验概率相等，均为1/3
library(MASS) # 加载R包
fit2 <- qda(type ~ Sepal.Length+Sepal.Width+Petal.Length+Petal.Width,data=na.omit(mydata),prior=c(1,1,1)/3) # prior用于指定先验概率，na.omit()去除缺失值

3. 结果可视化

你可以使用简单的plot()函数来可视化判别分析的结果（R语言入门之基本绘图），这里的横纵坐标分别代表前两个线性判别函数，每个观察点可以通过组别来区分。

# 使用前两个线性判别函数绘制散点图
plot(fit) # 用lda的拟合结果绘图

从上图我们可以看出0类和其它两类区分很明显，但是1和2的区分并不是非常完美，说明1和2这两类可能会有出错的可能。

接下来的代码根据第一线性判别函数绘制的每一类的柱状图和概率密度曲线：

# 就第一线性判别函数绘制各类的柱状图和概率密度曲线
plot(fit, dimen=1, type="both") # 使用lda的拟合结果

从上图我们可以看出，0组和其它两个组在X轴上没有重叠部分，但是1和2是由部分重叠的，说明0组和其它两组很容易分开，但是1和2组并不容易区分。这和散点图看出来的结果是一致的。

4. 假设条件的检验

在进行判别分析时，我们必须牢记：当数据是同方差时，我们可以使用线性判别函数；如果各组数据异方差，则使用二次判别函数更为准确。关于多元数据的异方差检验可以使用manova()函数。

关于判别分析的内容就讲到这里，咱们下期再见！