%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%应用多元统计分析 2010.04.02
%%%%Fisher判别函数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

> oil<-read.csv("Crudeoiltrans.csv")
> head(oil)
   X1       X2        X3        X4    X5   Group
1 3.9 7.141428 0.4472136 0.1416431 12.19 Wilhelm
2 2.7 7.000000 0.2645751 0.1400560 12.23 Wilhelm
3 2.8 6.000000 0.5477226 0.1428571 11.30 Wilhelm
4 3.1 6.708204 0.2828427 0.1388889 13.01 Wilhelm
5 3.5 6.782330 0.3162278 0.1280410 12.63 Wilhelm
6 3.9 6.557439 0.2645751 0.1600000 10.42 Wilhelm

> oilsweep<-sweep(oil[,1:5],2,apply(oil[1:5],2,mean)) #数据中心化
> e1<-matrix(c(.312,-.710,2.764,11.809,-.235),ncol=1) 
> y1<-as.matrix(oilsweep)%*%e1                #第一判别函数值
> e2<-matrix(c(.169,-.245,-2.046,-24.453,-.378),ncol=1)
> y2<-as.matrix(oilsweep)%*%e2                  #第二判别函数值
> oilfisher<-matrix(cbind(y1,y2,oil[,6]),ncol=3)

> m1<-apply(oilfisher[1:7,1:2],2,mean) 
> m2<-apply(oilfisher[8:18,1:2],2,mean)
> m3<-apply(oilfisher[19:56,1:2],2,mean)
> m<-matrix(rbind(m1,m2,m3),ncol=2)
> plot(c(-6,4), c(-3,3), xlab="Y1", ylab="Y2", type="n")
>  points(m,pch=19)
>  points(oilfisher[1:7,1:2],pch=21)
>  points(oilfisher[8:18,1:2],pch=22)
>  points(oilfisher[19:56,1:2],pch=24)