第七周-第八周-第九周
(1)示例2
am<-as.numeric(matrix(c(-1,-1),1))
bm<-as.numeric(matrix(c(2,1),1))
sp<-matrix(c(7.3,-1.1,-1.1,4.8),2,2)
(a=(am-bm)%*%solve(sp))
(y1=a%%am);(y2=a%%bm) #y1>y2
(y1=a%%am);(y2=a%%bm) #y1>y2
(y0=(y1+y2)/2)
(y0=as.vector(y0))
(x0=c(0,1))
(y=a%*%x0)
(y=as.vector(y))
y1>y2时，若y>y0,"判为第1类"
ifelse(y>y0,"判为第1类",
ifelse(y<y0,"判为第2类","无法判定"))
(2)
rm(list = ls())
a<-read.table("clipboard",header=T)#将数据mvexec4.xls中E6.3的数据处理并复制在剪贴板
attach(a)
线性判别分析
plot(x1,x2,pch=G,col=G)
text(x1,x2,G,adj=-0.5)
library(MASS)
(ld=lda(G~x1+x2))
z=predict(ld)
anewG=zclass
t<-table(G,a$newG)
sum(diag(prop.table(t)))#符合率0.8064516
距离分析
qd=qda(G~x1+x2)
z=predict(qd)
anewG=zclass
t<-table(G,a$newG)
sum(diag(prop.table(t)))#符合率.8387097
detach(a)
绘制线性判别函数图
(Y1=ldscaling[1]*ldmeans[1,1]+
ldscaling[2]*ldmeans[1,2])
第一类的两个指标的均值的“降维”(一维)
-11.87429
(Y2=ldscaling[1]*ldmeans[2,1]+
ldscaling[2]*ldmeans[2,2])
第二类的两个指标的均值的“降维”(一维)
-9.720089
Y1<Y2
(Y0=(Y1+Y2)/2)
abline(a=Y0/ldscaling[2], b=-ldscaling[1]/ld$scaling[2])
(3)
rm(list=ls())
b<-read.table("clipboard",header=T)#将数据mvexec4.xls中E6.4的数据忽略第一列并复制在剪贴板
attach(b)
线性判别分析
(ld=lda(G~x1+x2+x3+x4))
z=predict(ld)
bnewG=zclass
t<-table(G,b$newG)
sum(diag(prop.table(t)))#符合率 0.875
距离分析
(qd=qda(G~x1+x2+x3+x4))
z=predict(qd)
bnewG=zclass
t<-table(G,b$newG)
sum(diag(prop.table(t)))#符合率1
detach(b)
预测
p1<-predict(ld,data.frame(x1=78.3563,x2=0.8895,x3=1.8001,x4=14.1022))
p1$class
线性判别法判定其种类为ST公司，即在2009年该公司有可能会陷入财务困境
p2<-predict(qd,data.frame(x1=78.3563,x2=0.8895,x3=1.8001,x4=14.1022))
p2$class
距离判别法判定其种类为ST公司，即在2009年该公司有可能会陷入财务困境
第3题
方法一：使用MASS包的lda函数设置先验概率
rm(list=ls()) #清空内存
a<-read.table("clipboard",header=T)#将数据mvexec4.xls中E6.3的数据处理并复制在剪贴板
library(MASS)
先验概率相等的模型
(B1=lda(G~x1+x2,prior=c(1,1)/2))
attach(a)
P1=predict(B1)#判类
P1post#后验概率 round(P1post,3)#后验概率,显示小数点后3位
cbind(G,newG=P1class,后验概率=round(P1post,3))
矩阵的三列：真实类，判别类,后验概率
(t1=table(G,newG=P1$class))#混淆矩阵
(pr1=sum(diag(t1))/sum(t1))#判对率
先验概率“各总体样本数占样本总数的比例”的判别模型
(B2=lda(G~x1+x2,prior=c(15,16)/31))
P2=predict(B2)#判类
cbind(G,newG=P2class,后验概率=round(P2post,3))
矩阵的三列：真实类，判别类,后验概率
(t2=table(G,newG=P2$class))#混淆矩阵
(pr2=sum(diag(prop.table(t2))))#判对率0.8064516
detach(a)
方法二：使用WMDB包的dbayes函数设置先验概率
install.packages("WMDB")
library(WMDB)
dbayes(a[1:2],as.factor(G),p=rep(1/2,2))
先验概率相等(也是参数p的缺省值)的判别模型
dbayes(a[1:2],as.factor(G),p=c(15,16)/31)
先验概率不等的判别模型
第4题
方法一
rm(list=ls()) #清空内存
b<-read.table("clipboard",header=T)#将数据mvexec4.xls中E6.4的数据忽略第一列并复制在剪贴板
attach(b)
search()
(b1=lda(G~x1+x2+x3+x4,prior=c(7,8)/15))
方法二
dbayes(b[1:4],as.factor(G),p=rep(1/4,4))
鸢尾花习题
(d=iris)
attach(d)
线性判别
install.packages("MASS")
library(MASS)
ld1=lda(Species~.,data=d)
ld1lev attributes(ld1) pld=predict(ld1) #“预测”，即判别 data.frame(G=Species,newG=pldclass) #真实分类对比判别分类
(ta1=table(G=Species,newG=pld$class)) #混淆矩阵
addmargins(ta1) #添加行和、列和，总和
sum(diag(prop.table(ta1)))
判对率 0.98
二次判别
(qd=qda(Species~.,data=d))
Zqd=predict(qd) #“预测”，即判别
(ta3=table(G=Species,newG=Zqd$class)) #混淆矩阵
sum(diag(ta3/sum(ta3)))
判对率 0.98
先验概率相等的Bayes判别
(B=lda(Species~.,data=d,prior=c(1,1,1)/3))
ZB=predict(B) #“预测”，即判别
(ta2=table(G=Species,newG=Zld$class)) #混淆矩阵
sum(diag(ta2/sum(ta2)))
判对率 0.98
install.packages("WMDB")
library(WMDB)
dbayes函数判别
d[,6]=rep(1:3,each=50)
dbayes(d[,1:4],as.factor(d[,6]))
判对率 0.9666667
马氏距离判别
wmd(d[,1:4],as.factor(d[,6]))
