判别分析概念:
判别分析(Discriminat Analysis)是多元分析中用于判别样本所属类型的一种统计分析方法。
线性判别函数
小结:
1.判别分析方法是按一直所属组的样本确定判别函数,制定判别规则,然后再判断每一个新样品应属于哪一类。
2.常用的判别方法有Fisher判别、距离判别、贝叶斯判别等,每个方法根据其出发点不同各有其特点。
3.Fisher类判别对判别变量的分布类型并无要求,而Bayes类判别要变量的分布类型。因此,Fisher类判别较Bayes类判别简单一些。
4.当两个总体时,若它们的协方差矩阵相同,则距离判别和Fisher判别等价。当变量服从正态分布时,它们还和Bayes判别等价。
5.判别分析中的各种误判的后果允许看作时相同的,通常将犯第一类错误的后果看得更严重些。
library(openxlsx) #加载数据包
## Warning: package 'openxlsx' was built under R version 4.0.5
d6.3 = read.xlsx('mvstats5.xlsx','d6.3');d6.3
## Q C P G3
## 1 8.3 4.0 29 1
## 2 9.5 7.0 68 1
## 3 8.0 5.0 39 1
## 4 7.4 7.0 50 1
## 5 8.8 6.5 55 1
## 6 9.0 7.5 58 2
## 7 7.0 6.0 75 2
## 8 9.2 8.0 82 2
## 9 8.0 7.0 67 2
## 10 7.6 9.0 90 2
## 11 7.2 8.5 86 2
## 12 6.4 7.0 53 2
## 13 7.3 5.0 48 2
## 14 6.0 2.0 20 3
## 15 6.4 4.0 39 3
## 16 6.8 5.0 48 3
## 17 5.2 3.0 29 3
## 18 5.8 3.5 32 3
## 19 5.5 4.0 34 3
## 20 6.0 4.5 36 3
G3 = d6.3[,4]
Q = d6.3[,1]
C = d6.3[,2]
P = d6.3[,3] #绑定数据
library(MASS)
## Warning: package 'MASS' was built under R version 4.0.5
#(1).先验概率相等:q1=q2=q3=1/3,此时判别函数类似于Fisher线性判别函数
ld41 = lda(G3~Q+C+P,prior = c(1,1,1)/3) #先验概率相等的Bayes判别模型
#(2).先验概率不等:取q1=5/20,q2=8/20,q3=7/20,下面为先验概率不相等时的Eayes判别函数的系数
ld42=lda(G3~Q+C+P,prior=c(5,8,7)/20);ld42 #先验概率不相等的Bayes判别模型
## Call:
## lda(G3 ~ Q + C + P, prior = c(5, 8, 7)/20)
##
## Prior probabilities of groups:
## 1 2 3
## 0.25 0.40 0.35
##
## Group means:
## Q C P
## 1 8.400000 5.900000 48.200
## 2 7.712500 7.250000 69.875
## 3 5.957143 3.714286 34.000
##
## Coefficients of linear discriminants:
## LD1 LD2
## Q -0.81173396 0.88406311
## C -0.63090549 0.20134565
## P 0.01579385 -0.08775636
##
## Proportion of trace:
## LD1 LD2
## 0.7403 0.2597
#下面是两种结果比较:
Z1=predict(ld41);Z2=predict(ld42) #预测模型
data.frame(G3,ld41G=Z1$class,ld42G=Z2$class)#结果比较
## G3 ld41G ld42G
## 1 1 1 1
## 2 1 1 1
## 3 1 1 1
## 4 1 1 1
## 5 1 1 1
## 6 2 1 1
## 7 2 2 2
## 8 2 2 2
## 9 2 2 2
## 10 2 2 2
## 11 2 2 2
## 12 2 2 2
## 13 2 3 3
## 14 3 3 3
## 15 3 3 3
## 16 3 3 3
## 17 3 3 3
## 18 3 3 3
## 19 3 3 3
## 20 3 3 3
T1=table(G3,Z1$class);T1
##
## G3 1 2 3
## 1 5 0 0
## 2 1 6 1
## 3 0 0 7
T2=table(G3,Z2$class);T2
##
## G3 1 2 3
## 1 5 0 0
## 2 1 6 1
## 3 0 0 7
sum(diag(T1))/sum(T1) #diag对结果求对角线的数据和
## [1] 0.9
sum(diag(T1))/sum(T1)
## [1] 0.9
#ld41模型的后验概率
round(Z1$post*100,2)
## 1 2 3
## 1 98.26 0.56 1.19
## 2 79.42 20.57 0.01
## 3 93.72 4.31 1.97
## 4 65.37 33.71 0.91
## 5 90.52 9.44 0.05
## 6 92.78 7.21 0.00
## 7 0.33 86.32 13.34
## 8 17.75 82.25 0.01
## 9 18.47 81.05 0.48
## 10 0.28 99.70 0.02
## 11 0.22 99.69 0.09
## 12 11.12 77.98 10.90
## 13 29.18 32.50 38.32
## 14 0.08 0.02 99.90
## 15 1.21 2.27 96.52
## 16 7.94 24.27 67.79
## 17 0.01 0.04 99.95
## 18 0.14 0.28 99.58
## 19 0.10 0.43 99.47
## 20 1.38 2.52 96.10
#ld42模型的后验概率
round(Z2$post*100,2)
## 1 2 3
## 1 97.47 0.88 1.65
## 2 70.70 29.29 0.01
## 3 90.66 6.67 2.67
## 4 54.21 44.73 1.06
## 5 85.65 14.29 0.06
## 6 88.93 11.06 0.01
## 7 0.21 87.90 11.89
## 8 11.88 88.11 0.01
## 9 12.41 87.14 0.45
## 10 0.18 99.81 0.02
## 11 0.14 99.78 0.08
## 12 7.36 82.55 10.09
## 13 21.64 38.57 39.79
## 14 0.05 0.02 99.92
## 15 0.86 2.60 96.54
## 16 5.60 27.40 66.99
## 17 0.01 0.04 99.95
## 18 0.10 0.32 99.58
## 19 0.07 0.49 99.44
## 20 0.98 2.89 96.13
#后验概率给出了样品落在各个类的概率大小,这也是Bayes判别区别于Fisher判别的主要特点
predict(ld41,data.frame(Q=8,C=7.5,P=65))#ld41模型的判定
## $class
## [1] 2
## Levels: 1 2 3
##
## $posterior
## 1 2 3
## 1 0.300164 0.6980356 0.001800378
##
## $x
## LD1 LD2
## 1 -1.426693 -0.5046594
#根据建立的Bayes判别函数,带入预测数据,判断新样品属于第2类,即该产品实际上属于平销,但属于平销的概率仅69.8%。
predict(ld42,data.frame(Q=8,C=7.5,P=65))#ld42模型的判定,先验概率不相等的Bayes判别模型
## $class
## [1] 2
## Levels: 1 2 3
##
## $posterior
## 1 2 3
## 1 0.2114514 0.786773 0.001775594
##
## $x
## LD1 LD2
## 1 -1.537069 -0.1367865
#根据建立的Bayes判别函数,带入预测数据,判断新样品属于第2类,即该产品实际上属于平销,但属于平销的概率仅为78.7%,从这也可以看出考虑与不考虑先验概率对模型的判别效果还是有影响的。
##案例6企业财务状况的判别分析
Case6 = read.xlsx('mvcase5.xlsx','Case6');Case6
## G CF_TD NI_TA CA_CL CA_NS
## 1 1 0.51 0.10 2.49 0.54
## 2 1 0.08 0.02 2.01 0.53
## 3 1 0.38 0.11 3.27 0.35
## 4 1 0.19 0.05 2.25 0.33
## 5 1 0.32 0.07 4.24 0.63
## 6 1 0.31 0.05 4.45 0.69
## 7 1 0.12 0.05 2.52 0.69
## 8 1 -0.02 0.02 2.05 0.35
## 9 1 0.22 0.08 2.35 0.40
## 10 1 0.17 0.07 1.80 0.52
## 11 1 0.15 0.05 2.17 0.55
## 12 1 -0.10 -0.01 2.50 0.58
## 13 1 0.14 -0.03 0.46 0.26
## 14 1 0.14 0.07 2.61 0.52
## 15 1 0.15 0.06 2.23 0.56
## 16 1 0.16 0.05 2.31 0.20
## 17 1 0.29 0.06 1.84 0.38
## 18 1 0.54 0.11 2.33 0.48
## 19 1 -0.33 -0.09 3.01 0.47
## 20 1 0.48 0.09 1.24 0.18
## 21 1 0.56 0.11 4.29 0.44
## 22 1 0.20 0.08 1.99 0.30
## 23 1 0.47 0.14 2.92 0.45
## 24 1 0.17 0.04 2.45 0.14
## 25 1 0.58 0.04 5.06 0.13
## 26 2 -0.45 -0.41 1.09 0.45
## 27 2 -0.56 -0.31 1.51 0.16
## 28 2 0.06 0.02 1.01 0.40
## 29 2 -0.07 -0.09 1.45 0.26
## 30 2 -0.10 -0.09 1.56 0.67
## 31 2 -0.14 -0.07 0.71 0.28
## 32 2 0.04 0.01 1.50 0.71
## 33 2 -0.06 -0.06 1.37 0.40
## 34 2 0.07 -0.01 1.37 0.34
## 35 2 -0.13 -0.14 1.42 0.43
## 36 2 -0.23 -0.30 0.33 0.18
## 37 2 0.07 0.02 1.31 0.25
## 38 2 0.01 0.00 2.15 0.70
## 39 2 -0.28 -0.23 1.19 0.66
## 40 2 0.15 0.05 1.88 0.27
## 41 2 0.37 0.11 1.99 0.38
## 42 2 -0.08 -0.08 1.51 0.42
## 43 2 0.05 0.03 1.68 0.95
## 44 2 0.01 0.00 1.26 0.60
## 45 2 0.12 0.11 1.14 0.17
## 46 2 -0.28 -0.27 1.27 0.51
plot(Case6[,2:5],gap = 0)
library(MASS)
ld = lda(G~.,data = Case6);ld #线性判别
## Call:
## lda(G ~ ., data = Case6)
##
## Prior probabilities of groups:
## 1 2
## 0.5434783 0.4565217
##
## Group means:
## CF_TD NI_TA CA_CL CA_NS
## 1 0.23520000 0.05560000 2.593600 0.426800
## 2 -0.06809524 -0.08142857 1.366667 0.437619
##
## Coefficients of linear discriminants:
## LD1
## CF_TD -0.6291667
## NI_TA -4.4458516
## CA_CL -0.8892843
## CA_NS 1.1844801
plot(ld)
Zld = predict(ld) #线性判别模型预测模型
data.frame(Case6$G,Zld$class,round(Zld$x,3))
## Case6.G Zld.class LD1
## 1 1 1 -1.013
## 2 1 1 0.028
## 3 1 1 -1.895
## 4 1 1 -0.625
## 5 1 1 -2.210
## 6 1 1 -2.230
## 7 1 1 -0.395
## 8 1 1 -0.158
## 9 1 1 -0.783
## 10 1 1 -0.076
## 11 1 1 -0.268
## 12 1 1 -0.102
## 13 1 2 1.271
## 14 1 1 -0.778
## 15 1 1 -0.354
## 16 1 1 -0.813
## 17 1 1 -0.308
## 18 1 1 -1.005
## 19 1 1 -0.185
## 20 1 1 -0.265
## 21 1 1 -2.808
## 22 1 1 -0.569
## 23 1 1 -1.655
## 24 1 1 -0.971
## 25 1 1 -3.562
## 26 2 2 2.997
## 27 2 2 1.904
## 28 2 2 0.776
## 29 2 2 0.790
## 30 2 2 1.196
## 31 2 2 1.426
## 32 2 2 0.764
## 33 2 2 0.887
## 34 2 2 0.512
## 35 2 2 1.278
## 36 2 2 2.725
## 37 2 2 0.325
## 38 2 2 0.238
## 39 2 2 2.249
## 40 2 1 -0.342
## 41 2 1 -0.715
## 42 2 2 0.888
## 43 2 2 0.793
## 44 2 2 0.911
## 45 2 1 -0.050
## 46 2 2 2.178
tab1 = table(Case6$G,Zld$class);tab1
##
## 1 2
## 1 24 1
## 2 3 18
sum(diag(tab1))/sum(tab1) #准确率计算
## [1] 0.9130435
addmargins(tab1) #addmargins(table,margins) 将概述边margins放入表中
##
## 1 2 Sum
## 1 24 1 25
## 2 3 18 21
## Sum 27 19 46
qd = qda(G~.,data = Case6);qd #二次判别模型
## Call:
## qda(G ~ ., data = Case6)
##
## Prior probabilities of groups:
## 1 2
## 0.5434783 0.4565217
##
## Group means:
## CF_TD NI_TA CA_CL CA_NS
## 1 0.23520000 0.05560000 2.593600 0.426800
## 2 -0.06809524 -0.08142857 1.366667 0.437619
Zqd = predict(qd)
#data.frame(Case6$G,Zqd$class,round(Zqd$post,3)*100)
tab2 = table(Case6$G,Zqd$class);tab2
##
## 1 2
## 1 24 1
## 2 2 19
sum(diag(tab2))/sum(tab2)
## [1] 0.9347826
addmargins(tab2)
##
## 1 2 Sum
## 1 24 1 25
## 2 2 19 21
## Sum 26 20 46