（数据科学学习手札17）线性判别分析的原理简介&Python与R实现
之前数篇博客我们比较了几种具有代表性的聚类算法，但现实工作中，最多的问题是分类与定性预测，即通过基于已标注类型的数据的各显著特征值，通过大量样本训练出的模型，来对新出现的样本进行分类，这也是机器学习中最多的问题，而本文便要介绍分类算法中比较古老的线性判别分析：
线性判别
最早提出合理的判别分析法者是R.A.Fisher（1936），Fisher提出将线性判别函数用于花卉分类上，将花卉的各种特征利用线性组合方法变成单变量值，即将高维数据利用线性判别函数进行线性变化投影到一条直线上，再利用单值比较方法来对新样本进行分类，主要步骤如下：
Step1：求线性判别函数；
Step2：计算判别界值；
Step3：建立判别标准（这里与模糊分类中的隶属度有些相似，即离哪一类的投影中心最近，就将样本判别为哪一类）
下面分别利用Python，R，基于著名的花卉分类数据集iris进行演示：
Python
我们利用sklearn包中封装的LinearDiscriminantAnalysis对iris构建线性判别模型，因为LDA实际上是将高维数据尽可能分开的投影到一条直线上，因此LDA也可以对特定数据进行降维转换：
'''Fisher线性判别分析''' import numpy as np import matplotlib.pyplot as plt from sklearn import datasets from sklearn.discriminant_analysis import LinearDiscriminantAnalysis from matplotlib.pyplot import style from sklearn.model_selection import train_test_split style.use('ggplot') iris = datasets.load_iris() X = iris.data y = iris.target '''展示LDA的降维功能''' target_names = iris.target_names '''设置压缩到1维''' lda = LinearDiscriminantAnalysis(n_components=1) '''利用线性判别函数将四维的样本数据压缩到一条直线上''' X_r2 = lda.fit(X,y).transform(X) X_Zero = np.zeros(X_r2.shape) '''绘制降维效果图''' for c,i,target_names in zip('ryb',[0,1,2],target_names): plt.scatter(X_r2[y == i],X_Zero[y == i],c=c,label=target_names,s=5) plt.legend() plt.grid()
降维后的效果图如下：
下面正式对iris数据集进行LDA分类，这里用到一个新的方法，是sklearn.model_selection.train_test_split，它的作用是根据设置的训练集与测试集的比例进行随机分割，我们利用从样本集中分割的7成数据作为训练集，3成数据进行测试，过程及结果如下：
'''利用sklearn自带的样本集划分方法进行分类，这里选择训练集测试集73开''' X_train,X_test,y_train,y_test = train_test_split(X,y,test_size=0.3) '''搭建LDA模型''' lda = LinearDiscriminantAnalysis(n_components=1) '''利用分割好的训练集进行模型训练并对测试集进行预测''' ld = lda.fit(X_train,y_train).predict(X_test) '''比较预测结果与真实分类结果''' print(np.array([ld,y_test])) '''打印正确率''' print('正确率：',str(lda.score(X_test,y_test)))
结果如下：
可以看出，在iris上取得了非常高的准确率。
R
在R中做LDA需要用到MASS包中的lda(formula~feature1+feature2+...+featuren,data=df)，其中formula表示数据集中表示分类标注的列，右边的各种feature表示将要使用到的分类特征值，也即是构建线性判别函数要用到的基础变量，data指保存全部数据的数据框，具体过程如下：
> #Fisher线性判别 > rm(list=ls()) > library(MASS) > data(iris) > data <- iris > data$Species = as.character(data$Species) > #创造类别变量 > data$type[data$Species == 'setosa'] = 1 > data$type[data$Species == 'versicolor'] = 2 > data$type[data$Species == 'virginica'] = 3 > #利用简单随机抽样将样本集划分为训练集与验证集 > sam <- sample(1:length(data[,1]),105) > train <- data[sam,] > test <- data[-sam,] > #根据样本数据创建线性判别分析模型 > ld <- lda(type~Sepal.Length+Sepal.Width+Petal.Length+Petal.Width,data=train) > #将样本集作为验证集求分类结果 > Z <- predict(ld) > #保存预测类别 > newType <- Z$class > #与真实分类结果进行比较 > cbind(train$type,Z$x,newType) LD1 LD2 newType 74 2 -2.0419740 -1.16551052 2 93 2 -1.2286978 -1.34868982 2 124 3 -4.2056288 -0.39378496 3 43 1 7.0499296 -0.38479505 1 141 3 -6.3297144 1.43240686 3 112 3 -5.2009805 -0.39665336 3 75 2 -1.2076317 -0.53028083 2 83 2 -0.8607055 -1.04255122 2 128 3 -3.7958024 0.38283281 3 134 3 -3.5521695 -0.79338125 2 84 2 -4.1930475 -0.86304028 3 80 2 0.2455854 -1.48598565 2 16 1 8.9754828 2.99336222 1 35 1 6.7862663 -0.74053169 1 64 2 -2.4399568 -0.53150309 2 79 2 -2.4024680 -0.24830622 2 142 3 -4.9647448 1.48158101 3 54 2 -2.2217845 -1.93554530 2 14 1 7.3092968 -0.99420385 1 119 3 -8.8253015 -0.67042900 3 71 2 -3.3701925 0.94839271 3 44 1 6.2098650 0.99672773 1 29 1 7.6841973 0.08156876 1 62 2 -1.7099434 0.15332464 2 59 2 -1.6648086 -0.67372500 2 77 2 -2.3932456 -0.84052244 2 27 1 6.6078477 0.36272030 1 101 3 -7.2690741 1.94879260 3 91 2 -2.1739793 -1.53761667 2 40 1 7.4330661 0.03443993 1 103 3 -5.9833783 0.46642759 3 121 3 -5.9399677 1.45493705 3 51 2 -1.3784535 0.24058778 2 88 2 -2.5153892 -2.12899905 2 6 1 7.5334551 1.60652494 1 148 3 -4.7085798 0.61397476 3 76 2 -1.3919721 -0.13234144 2 10 1 7.0733099 -0.92823629 1 7 1 7.0416936 0.27174257 1 104 3 -5.1966351 -0.20937105 3 24 1 6.0120553 0.24387476 1 137 3 -6.0383747 2.20984770 3 61 2 -1.2257443 -3.03469429 2 41 1 7.6466593 0.57623503 1 26 1 6.4775174 -1.04708182 1 70 2 -1.0443574 -1.74662921 2 31 1 6.5351351 -0.78766052 1 85 2 -2.5818237 0.01276122 2 36 1 7.5972776 -0.33972390 1 4 1 6.6085362 -0.73929708 1 130 3 -4.2970292 -0.42496657 3 95 2 -1.8418994 -0.99664464 2 111 3 -4.1644831 1.17871159 3 117 3 -4.7101562 0.09594447 3 48 1 6.9765285 -0.43315849 1 106 3 -7.0303766 0.13158736 3 69 2 -3.5167119 -2.05931691 2 131 3 -5.9675456 -0.52249341 3 92 2 -2.0719645 -0.22536449 2 22 1 7.3872921 1.18564382 1 39 1 6.6977209 -0.90199151 1 125 3 -5.3082628 1.33894916 3 47 1 7.9455957 1.02129249 1 108 3 -5.9474167 -0.54626897 3 123 3 -7.2281863 -0.62126561 3 25 1 6.4877846 -0.15448692 1 96 2 -0.9672993 -0.40896608 2 149 3 -5.4267988 2.11763536 3 58 2 -0.1967945 -1.90480909 2 100 2 -1.4146638 -0.69091758 2 107 3 -4.3326495 -0.90276305 3 89 2 -1.1216984 -0.17330958 2 67 2 -2.4633370 0.01193815 2 45 1 6.7958450 1.25408059 1 138 3 -4.5932952 0.35495424 3 30 1 6.6519961 -0.52865075 1 17 1 8.3010066 1.79668640 1 9 1 6.3297286 -1.20813011 1 135 3 -4.6952607 -1.73516107 3 18 1 7.5140148 0.52828312 1 145 3 -6.4564371 2.08976755 3 68 2 -0.6703941 -1.51304115 2 12 1 7.0634482 -0.01186583 1 23 1 8.4484974 0.79139590 1 20 1 7.8504400 1.25653745 1 15 1 9.4800594 1.72576966 1 147 3 -5.0367690 -0.77081719 3 46 1 6.4557628 -0.76347342 1 60 2 -1.7902528 -0.66467281 2 19 1 7.8221244 1.15898751 1 144 3 -6.3829867 1.36026786 3 1 1 7.8010583 0.34057853 1 90 2 -1.8695758 -1.41834884 2 53 2 -2.2846205 0.07502495 2 115 3 -6.4317784 0.89801781 3 146 3 -5.4512238 1.17626549 3 78 2 -3.3451866 0.14511863 2 99 2 0.3864166 -1.31670824 2 118 3 -6.0412297 2.34012590 3 38 1 8.1457196 0.41229523 1 139 3 -3.6631579 0.43078471 3 105 3 -6.4339941 0.70414177 3 132 3 -4.7729922 2.10651473 3 110 3 -6.3994587 2.67334311 3 72 2 -0.9858025 -0.64502336 2 > #打印混淆矩阵 > (tab <- table(newType,train$type)) newType 1 2 3 1 35 0 0 2 0 33 1 3 0 2 34 > #显示正确率 > cat('Accuracy:',sum(diag(tab))/length(train[,1])) Accuracy: 0.9714286
> #将验证集代入训练好的模型中计算近似泛化误差 > T <-predict(ld,test) > #与真实分类结果进行比较 > cbind(test$type,T$x,T$class) LD1 LD2 2 1 6.8020498 -0.95158956 1 3 1 7.2276597 -0.38602966 1 5 1 7.9179194 0.59958829 1 8 1 7.3738228 0.03485147 1 11 1 8.1391093 0.80900002 1 13 1 7.0298500 -1.13888262 1 21 1 7.2270204 -0.06187541 1 28 1 7.6684138 0.29262662 1 32 1 7.0367090 0.40861452 1 33 1 9.0120836 1.65651141 1 34 1 9.2707624 2.14912000 1 37 1 8.2299196 0.38647275 1 42 1 5.2371900 -2.52488607 1 49 1 8.0798659 0.80941155 1 50 1 7.3896062 -0.17620640 1 52 2 -1.6371815 0.52584233 2 55 2 -2.4742435 -0.55650250 2 56 2 -2.1822152 -0.88107904 2 57 2 -2.1911397 0.87747596 2 63 2 -1.2405414 -2.75931501 2 65 2 -0.3383635 -0.19420599 2 66 2 -1.1566243 0.12584525 2 73 2 -3.6967069 -1.47409522 2 81 2 -1.0878173 -1.95727554 2 82 2 -0.6088858 -2.09743977 2 86 2 -1.8089897 1.23238953 2 87 2 -2.0193315 0.17092876 2 94 2 -0.3136555 -2.16381886 2 97 2 -1.4304473 -0.47985972 2 98 2 -1.3261184 -0.52945776 2 102 3 -5.1726649 -0.29910341 3 109 3 -6.0478550 -1.34049085 3 113 3 -5.3935569 0.65782366 3 114 3 -5.6792727 -0.58064337 3 116 3 -5.4686329 1.64715619 3 120 3 -4.5946379 -2.29619567 3 122 3 -5.0183151 0.24310322 3 126 3 -4.9026833 0.37255836 3 127 3 -3.8968800 -0.08723483 3 129 3 -6.1746269 0.09473298 3 133 3 -6.4616705 0.28243757 3 136 3 -6.5857811 0.74428684 3 140 3 -4.9663213 0.96355072 3 143 3 -5.1726649 -0.29910341 3 150 3 -4.2980648 0.28857515 3 > #打印混淆矩阵 > (tab <- table(T$class,test$type)) 1 2 3 1 15 0 0 2 0 15 0 3 0 0 15 > #显示正确率 > cat('Accuracy:',sum(diag(tab))/length(test[,1])) Accuracy: 1 > #Fisher线性判别 > rm(list=ls()) > library(MASS) > data(iris) > data <- iris > data$Species = as.character(data$Species) > #创造类别变量 > data$type[data$Species == 'setosa'] = 1 > data$type[data$Species == 'versicolor'] = 2 > data$type[data$Species == 'virginica'] = 3 > #利用简单随机抽样将样本集划分为训练集与验证集 > sam <- sample(1:length(data[,1]),105) > train <- data[sam,] > test <- data[-sam,] > #根据样本数据创建线性判别分析模型 > ld <- lda(type~Sepal.Length+Sepal.Width+Petal.Length+Petal.Width,data=train) > #将样本集作为验证集求分类结果 > Z <- predict(ld) > #保存预测类别 > newType <- Z$class > #与真实分类结果进行比较 > cbind(train$type,Z$x,newType) LD1 LD2 newType 44 1 6.5356943 1.380095755 1 150 3 -5.0777607 0.353229454 3 48 1 7.2092979 -0.233679182 1 17 1 8.7810490 1.832128789 1 11 1 8.5106371 0.591364581 1 13 1 7.4261015 -0.949875862 1 82 2 -0.7168056 -1.805629506 2 19 1 8.2170070 0.870266056 1 123 3 -7.9020752 -1.233817765 3 23 1 8.7554448 0.980676590 1 59 2 -1.8788410 -0.915972039 2 70 2 -1.2367635 -1.516788339 2 146 3 -5.6852127 1.788579457 3 142 3 -5.0723482 2.082228707 3 115 3 -6.9597642 1.869494439 3 7 1 7.2188221 0.472577805 1 112 3 -5.6667730 -0.156736416 3 117 3 -5.4232144 -0.112017744 3 26 1 6.8888800 -0.860736454 1 111 3 -4.6405840 1.346460523 3 29 1 8.1491362 0.125965565 1 81 2 -1.2242456 -1.597693375 2 14 1 7.6286098 -0.598455540 1 110 3 -7.1444731 2.585994755 3 40 1 7.7798631 0.022078382 1 39 1 6.9782244 -0.504346286 1 68 2 -1.0104685 -1.671771426 2 20 1 8.0502422 1.117799195 1 31 1 6.7929797 -0.656064848 1 134 3 -4.1319678 -1.006493240 2 73 2 -4.0110318 -1.372974920 2 4 1 6.8227612 -0.537268125 1 108 3 -6.6972047 -1.217364668 3 54 2 -2.3695531 -1.376553027 2 148 3 -5.1575482 0.848139667 3 141 3 -6.8305530 1.864676411 3 109 3 -6.6173517 -1.422188092 3 50 1 7.7923809 -0.058826655 1 91 2 -2.7513769 -1.544578568 2 97 2 -1.8728790 -0.435815301 2 96 2 -1.4911207 -0.557876549 2 77 2 -2.5464899 -1.021564544 2 66 2 -1.2439344 0.003029162 2 145 3 -7.0771380 2.462076368 3 38 1 8.3218303 0.213545770 1 45 1 6.7744745 0.999164236 1 58 2 -0.3713519 -1.340382308 2 35 1 7.1622529 -0.552177665 1 36 1 8.1741719 -0.035844508 1 57 2 -2.7067243 0.719387605 2 27 1 6.9079284 0.551777520 1 132 3 -5.6031449 1.030132762 3 128 3 -4.3392145 0.561003821 3 87 2 -2.2635929 -0.006748767 2 88 2 -2.4886790 -1.851739918 2 118 3 -7.1004619 1.347087687 3 119 3 -9.4291409 -0.890616169 3 46 1 6.9234400 -0.316289540 1 33 1 9.0573501 1.063438766 1 10 1 7.4135836 -0.868970825 1 78 2 -3.6650096 0.105534552 2 32 1 7.6166932 0.640755164 1 139 3 -4.1962691 0.674830697 3 103 3 -6.5226603 0.373114540 3 43 1 7.2390794 -0.114882460 1 121 3 -6.4785762 1.623818440 3 15 1 10.1229028 1.482252025 1 125 3 -6.0718150 1.194903725 3 92 2 -2.5655637 -0.379597733 2 21 1 7.6071362 -0.210545218 1 116 3 -6.0199587 2.084095792 3 95 2 -2.2468978 -0.820309281 2 93 2 -1.3874483 -1.124059988 2 135 3 -5.6483660 -2.247095684 3 52 2 -1.9604386 0.420606745 2 37 1 8.8751646 0.414644969 1 127 3 -4.2307963 0.275427178 3 137 3 -6.8919261 2.468751536 3 98 2 -1.5631688 -0.569521563 2 72 2 -1.0384324 -0.432712540 2 80 2 0.2825956 -1.208391313 2 67 2 -3.1266046 0.070901690 2 120 3 -4.9979151 -2.051118149 3 114 3 -6.2027782 0.131952933 3 69 2 -3.4910411 -1.516772693 2 130 3 -4.8967334 -1.106964081 3 42 1 5.9270652 -1.555646363 1 64 2 -2.9521004 -0.683186676 2 86 2 -2.4035699 1.146743118 2 75 2 -1.3368410 -0.579461256 2 143 3 -5.8335378 0.090796722 3 140 3 -5.3380144 1.122071295 3 1 1 8.1663998 0.325667325 1 83 2 -1.0009115 -0.820471045 2 28 1 8.0234544 0.211840448 1 85 2 -3.3529324 0.080841383 2 89 2 -1.5995061 -0.127256512 2 18 1 7.9150690 0.642460485 1 147 3 -5.2724641 -0.214659306 3 105 3 -7.1968652 0.828583809 3 41 1 8.0580144 0.756287362 1 63 2 -1.1801689 -2.546513654 2 122 3 -5.6685504 0.829975549 3 107 3 -5.0735507 -0.234382628 3 47 1 8.0454637 0.692149004 1 > #打印混淆矩阵 > (tab <- table(newType,train$type)) newType 1 2 3 1 36 0 0 2 0 33 1 3 0 0 35 > #显示正确率 > cat('Accuracy:',sum(diag(tab))/length(train[,1])) Accuracy: 0.9904762> #将验证集代入训练好的模型中计算近似泛化误差 > T <-predict(ld,test) > #与真实分类结果进行比较 > cbind(test$type,T$x,T$class) LD1 LD2 2 1 7.2879346 -0.63805255 1 3 1 7.5785711 -0.12979200 1 5 1 8.1836634 0.52536908 1 6 1 7.7566121 1.39670067 1 8 1 7.6666992 0.02704823 1 9 1 6.5916877 -0.80793523 1 12 1 7.1842622 -0.07186911 1 16 1 9.2604596 2.57316476 1 22 1 7.6684840 1.23986044 1 24 1 6.3832248 0.56001189 1 25 1 6.4159345 -0.39844020 1 30 1 6.8102433 -0.45636309 1 34 1 9.5320477 1.66891133 1 49 1 8.3974733 0.59633443 1 51 2 -1.5423430 -0.14371955 2 53 2 -2.5494836 -0.23440252 2 55 2 -2.6250939 -0.47214778 2 56 2 -2.7716342 -0.95711830 2 60 2 -2.1825564 -0.15706564 2 61 2 -1.2921165 -2.34199387 2 62 2 -2.0187853 0.38256324 2 65 2 -0.4493874 0.22229672 2 71 2 -4.0485780 1.06926437 3 74 2 -2.5798663 -1.51150491 2 76 2 -1.4875258 -0.18673290 2 79 2 -2.8043766 -0.14370961 2 84 2 -4.8532177 -0.86952245 3 90 2 -2.1086982 -0.98708920 2 94 2 -0.3886155 -1.54008407 2 99 2 0.5024002 -0.51222584 2 100 2 -1.7471972 -0.52169018 2 101 3 -8.2981213 2.15538466 3 102 3 -5.8335378 0.09079672 3 104 3 -6.0360789 -0.40566699 3 106 3 -7.7496056 -0.41373390 3 113 3 -5.8377150 0.82345220 3 124 3 -4.5041691 -0.03313161 3 126 3 -5.6507584 -0.30162799 3 129 3 -6.8073348 0.34501073 3 131 3 -6.4535806 -0.88255921 3 133 3 -7.0586655 0.66180389 3 136 3 -6.8585571 0.75916772 3 138 3 -5.4059508 0.08768402 3 144 3 -7.1039586 1.41107423 3 149 3 -6.2415407 2.37464228 3 > #打印混淆矩阵 > (tab <- table(T$class,test$type)) 1 2 3 1 14 0 0 2 0 15 0 3 0 2 14 > #显示正确率 > cat('Accuracy:',sum(diag(tab))/length(test[,1])) Accuracy: 0.9555556
可以看出，和Python中的效果相差无几。
以上就是关于线性判别的基本内容，如有意见望提出。
