#title 판별분석(discriminant analysis) [[TableOfContents]] ==== 예제 데이터 ==== [http://www.kyobobook.co.kr/product/detailViewKor.laf?ejkGb=KOR&mallGb=KOR&barcode=9788983257000&orderClick=LAG&Kc=SETLBkserp11_15 EXCEL에 의한 조사방법 및 통계분석]를 이용함. {{{ cname <- c("ID", "구매브랜드", "연령","세대연수입", "세대사람수", "방문빈도", "거주년수") x = read.table("c:\\data\\disc.txt", col.names = cname) head(x) }}} attachment:판별분석/disc.txt {{{ > head(x) ID 구매브랜드 연령 세대연수입 세대사람수 방문빈도 거주년수 1 1 A 48 9000 4 5 6 2 2 A 58 8000 6 4 20 3 3 A 52 7000 6 4 12 4 4 A 63 7000 6 4 15 5 5 A 59 8000 4 6 6 6 6 A 38 11000 5 4 10 > }}} 구매브랜드 A, B별로 평균(mean)과 분산(var)을 알고 싶으면 각각 서브셋을 만든 후에 계산한다. {{{ x1 <- subset(x, 구매브랜드=="A") x2 <- subset(x, 구매브랜드=="B") }}} ==== 피셔의 판별함수 ==== {{{ library("MASS") rs<-lda(구매브랜드~연령+세대연수입+세대사람수+방문빈도+거주년수, data=x) rs }}} 결과 {{{ > rs Call: lda(구매브랜드 ~ 연령 + 세대연수입 + 세대사람수 + 방문빈도 + 거주년수, data = x) Prior probabilities of groups: A B 0.5 0.5 Group means: 연령 세대연수입 세대사람수 방문빈도 거주년수 A 55.0 8400 4.8 4.6 9.3 B 32.9 5400 2.3 1.9 2.5 Coefficients of linear discriminants: LD1 연령 -0.0982250865 세대연수입 -0.0008926856 세대사람수 -0.7287341389 방문빈도 -0.6346846738 거주년수 -0.0578973013 > }}} ==== 결과해석 ==== {{{ p<-predict(rs) x <- cbind(x, p$x) x }}} 결과 {{{ > x ID 구매브랜드 연령 세대연수입 세대사람수 방문빈도 거주년수 LD1 1 1 A 48 9000 4 5 6 -3.716870 2 2 A 58 8000 6 4 20 -5.439781 3 3 A 52 7000 6 4 12 -3.494566 4 4 A 63 7000 6 4 15 -4.748734 5 5 A 59 8000 4 6 6 -4.539345 6 6 A 38 11000 5 4 10 -4.845629 7 7 A 49 9000 5 5 9 -4.717521 8 8 A 62 8000 3 7 5 -4.682073 9 9 A 55 9000 5 3 5 -3.805913 10 10 A 66 8000 4 4 5 -3.899654 11 11 B 28 6000 2 2 2 4.518800 12 12 B 33 4000 3 3 3 4.391730 13 13 B 26 5000 3 1 2 5.513887 14 14 B 40 7000 2 2 1 2.505311 15 15 B 31 5000 2 1 3 5.693598 16 16 B 24 6000 2 1 3 5.488488 17 17 B 47 6000 3 1 1 2.616372 18 18 B 30 7000 1 2 3 4.100501 19 19 B 40 4000 1 2 2 5.854205 20 20 B 30 4000 4 4 5 3.207192 > }}} LD1값이 A브랜드를 구매한 사람은 음수이고, B브랜드를 구매한 사람은 양수다. 잘 판별되었다고 보면 된다. 다음을 실행해보면 계산된 결과가 15.46785를 기준으로 나뉜다는 것을 알 수 있다. {{{ calc <- with(x, 연령 * -0.0982250865 + 세대연수입 * -0.0008926856 + 세대사람수 * -0.7287341389 + 방문빈도 * -0.6346846738 + 거주년수 * -0.0578973013) x <- cbind(x, calc) calc2 <- with(x, LD1 - calc) x <- cbind(x, calc2) x }}} 결과 {{{ > x ID 구매브랜드 연령 세대연수입 세대사람수 방문빈도 거주년수 LD1 calc calc2 1 1 A 48 9000 4 5 6 -3.716870 -19.184718 15.46785 2 2 A 58 8000 6 4 20 -5.439781 -20.907629 15.46785 3 3 A 52 7000 6 4 12 -3.494566 -18.962415 15.46785 4 4 A 63 7000 6 4 15 -4.748734 -20.216583 15.46785 5 5 A 59 8000 4 6 6 -4.539345 -20.007193 15.46785 6 6 A 38 11000 5 4 10 -4.845629 -20.313477 15.46785 7 7 A 49 9000 5 5 9 -4.717521 -20.185369 15.46785 8 8 A 62 8000 3 7 5 -4.682073 -20.149922 15.46785 9 9 A 55 9000 5 3 5 -3.805913 -19.273761 15.46785 10 10 A 66 8000 4 4 5 -3.899654 -19.367502 15.46785 11 11 B 28 6000 2 2 2 4.518800 -10.949048 15.46785 12 12 B 33 4000 3 3 3 4.391730 -11.076119 15.46785 13 13 B 26 5000 3 1 2 5.513887 -9.953962 15.46785 14 14 B 40 7000 2 2 1 2.505311 -12.962538 15.46785 15 15 B 31 5000 2 1 3 5.693598 -9.774251 15.46785 16 16 B 24 6000 2 1 3 5.488488 -9.979361 15.46785 17 17 B 47 6000 3 1 1 2.616372 -12.851477 15.46785 18 18 B 30 7000 1 2 3 4.100501 -11.367347 15.46785 19 19 B 40 4000 1 2 2 5.854205 -9.613644 15.46785 20 20 B 30 4000 4 4 5 3.207192 -12.260657 15.46785 }}} 세대사람수로 판별할 수도 있다. {{{ cname <- c("ID", "구매브랜드", "연령","세대연수입", "세대사람수", "방문빈도", "거주년수") x = read.table("c:\\data\\disc.txt", col.names = cname) rs<-lda(세대사람수~연령+세대연수입+구매브랜드+방문빈도+거주년수, data=x) rs p<-predict(rs) x <- cbind(x, p$x) x }}} 결과 {{{ > rs Call: lda(세대사람수 ~ 연령 + 세대연수입 + 구매브랜드 + 방문빈도 + 거주년수, data = x) Prior probabilities of groups: 1 2 3 4 5 6 0.10 0.20 0.20 0.20 0.15 0.15 Group means: 연령 세대연수입 구매브랜드B 방문빈도 거주년수 1 35.00000 5500.000 1.00 2.00 2.50000 2 30.75000 6000.000 1.00 1.50 2.25000 3 42.00000 5750.000 0.75 3.00 2.75000 4 50.75000 7250.000 0.25 4.75 5.50000 5 47.33333 9666.667 0.00 4.00 8.00000 6 57.66667 7333.333 0.00 4.00 15.66667 Coefficients of linear discriminants: LD1 LD2 LD3 LD4 연령 0.0108341882 -0.0842979470 -0.0182453017 0.0079837572 세대연수입 -0.0005159674 0.0003233735 0.0005008073 0.0005652328 구매브랜드B -3.5186494778 -3.6997526685 1.3377932310 5.0233382564 방문빈도 -0.5541865385 -0.0619435929 -0.5560926546 0.7505362060 거주년수 0.5078180145 -0.1710574267 0.1361173174 0.1627906244 LD5 연령 -0.1117330322 세대연수입 -0.0005554468 구매브랜드B -3.7199816354 방문빈도 0.0469918970 거주년수 -0.0118410355 Proportion of trace: LD1 LD2 LD3 LD4 LD5 0.8430 0.0983 0.0578 0.0006 0.0004 > x ID 구매브랜드 연령 세대연수입 세대사람수 방문빈도 거주년수 LD1 LD2 LD3 LD4 LD5 1 1 A 48 9000 4 5 6 -0.1993729 2.06204706 -0.6506451 0.03737132 0.3220854 2 2 A 58 8000 6 4 20 8.0885751 -1.43716633 1.1278296 1.08050866 -0.4525645 3 3 A 52 7000 6 4 12 4.4769932 0.11370722 -0.3524444 -0.83495164 0.8680088 4 4 A 63 7000 6 4 15 6.1196233 -1.32674248 -0.1447908 -0.25875844 -0.3965776 5 5 A 59 8000 4 6 6 -0.1184160 0.74945251 -1.9082434 0.31049609 -0.3045392 6 6 A 38 11000 5 4 10 1.2458091 2.92948751 1.6339845 0.98862558 0.2341661 7 7 A 49 9000 5 5 9 1.3349153 1.46457684 -0.2605384 0.53372695 0.1748293 8 8 A 62 8000 3 7 5 -1.1479180 0.60567250 -2.6551893 0.92219294 -0.5809054 9 9 A 55 9000 5 3 5 0.4770215 1.76690605 0.1977058 -1.57060541 -0.5421886 10 10 A 66 8000 4 4 5 0.5579784 0.45431149 -1.0598925 -1.29748065 -1.1688132 11 11 B 28 6000 2 2 2 -2.7555165 -0.05180681 0.6734409 0.30256502 0.4094934 12 12 B 33 4000 3 3 3 -1.7157794 -1.35304466 -0.8393757 0.12534510 0.9968727 13 13 B 26 5000 3 1 2 -1.7070310 -0.14464087 0.7652168 -1.02917147 1.1414143 14 14 B 40 7000 2 2 1 -3.6492916 -0.56895121 0.8191873 0.80081225 -1.4749088 15 15 B 31 5000 2 1 3 -1.1450420 -0.73718803 0.8101076 -0.82646206 0.5709081 16 16 B 24 6000 2 1 3 -1.7368487 0.17627114 1.4386320 -0.31711559 0.7975926 17 17 B 47 6000 3 1 1 -2.5032984 -1.42046679 0.7467555 -0.45907043 -1.7485851 18 18 B 30 7000 1 2 3 -2.7419975 -0.06808659 1.2738749 1.04655592 -0.3812606 19 19 B 40 4000 1 2 2 -1.5935715 -1.71012927 -0.5471174 -0.73209543 0.1795906 20 20 B 30 4000 4 4 5 -1.2868324 -1.50420926 -1.0684978 1.17751128 1.3553816 }}} ==== iris ==== {{{ tmp <- iris library("MASS") rs<-lda(Species~Sepal.Length+Sepal.Width+Petal.Length+Petal.Width, data=tmp) rs }}} 결과 {{{ > rs Call: lda(Species ~ Sepal.Length + Sepal.Width + Petal.Length + Petal.Width, data = tmp) Prior probabilities of groups: setosa versicolor virginica 0.3333333 0.3333333 0.3333333 Group means: Sepal.Length Sepal.Width Petal.Length Petal.Width setosa 5.006 3.428 1.462 0.246 versicolor 5.936 2.770 4.260 1.326 virginica 6.588 2.974 5.552 2.026 Coefficients of linear discriminants: LD1 LD2 Sepal.Length 0.8293776 0.02410215 Sepal.Width 1.5344731 2.16452123 Petal.Length -2.2012117 -0.93192121 Petal.Width -2.8104603 2.83918785 Proportion of trace: LD1 LD2 0.9912 0.0088 > }}} Proportion of trace가 LD1이 0.9912임. LD1만 써도 되것다. {{{ p <- predict(rs) tmp <- cbind(tmp, p$x) tmp[5:7] plot(rs) plot(rs, dimen=1, type="both") }}} attachment:판별분석/lda01.png x축이 LD1이고, y축이 LD2다. 뭉쳐있는 부분을 시각적으로 확인 할 수 있다. attachment:판별분석/lda02.png 히스토그림으로 그룹별로 데이터의 분포가 어떠한지 알 수 있다. {{{ mean(tmp$LD1[tmp$Species == "setosa"]) mean(tmp$LD1[tmp$Species == "versicolor"]) mean(tmp$LD1[tmp$Species == "virginica"]) }}} 결과 {{{ > mean(tmp$LD1[tmp$Species == "setosa"]) [1] 7.6076 > mean(tmp$LD1[tmp$Species == "versicolor"]) [1] -1.825049 > mean(tmp$LD1[tmp$Species == "virginica"]) [1] -5.78255 }}} ==== 참고자료 ==== * http://www.statmethods.net/advstats/discriminant.html * attachment:판별분석/discrim.pdf