为了研究世界各国森林、草原资源的分布规律,共抽取了21个国家的数据,每份数据包括4项指标,原始数据见表5-2。试用该数据对国别进行聚类分析。[大谦MATLAB,dqmatlab点com]
表5-2 原始数据表
| 国 别 | 森林面积(万公顷) | 森林覆盖率(%) | 林木蓄积量(亿立方米) | 草原面积(万公顷) |
|---|---|---|---|---|
| 中 国 | 11978 | 12.5 | 93.5 | 31908 |
| 美 国 | 28446 | 30.4 | 202.0 | 23754 |
| 日 本 | 2501 | 67.2 | 24.8 | 58 |
| 德 国 | 1028 | 28.4 | 14.0 | 599 |
| 英 国 | 210 | 8.6 | 1.5 | 1147 |
| 法 国 | 1458 | 26.7 | 16.0 | 1288 |
| 意 大 利 | 635 | 21.1 | 3.6 | 514 |
| 加 拿 大 | 32613 | 32.7 | 192.8 | 2385 |
| 澳大利亚 | 10700 | 13.9 | 10.5 | 45190 |
| 俄 罗 斯 | 92000 | 41.1 | 841.5 | 37370 |
| 捷 克 | 458 | 35.8 | 8.9 | 168 |
| 波 兰 | 868 | 27.8 | 11.4 | 405 |
| 匈 牙 利 | 161 | 17.4 | 2.5 | 129 |
| 克罗地亚 | 929 | 36.3 | 11.4 | 640 |
| 罗马尼亚 | 634 | 26.7 | 11.3 | 447 |
| 保加利亚 | 385 | 34.7 | 2.5 | 200 |
| 印 度 | 6748 | 20.5 | 29.0 | 1200 |
| 印度尼西亚 | 2180 | 84.0 | 33.7 | 1200 |
| 尼日利亚 | 1490 | 16.1 | 0.8 | 2090 |
| 墨 西 哥 | 4850 | 24.6 | 32.6 | 7450 |
| 巴 西 | 57500 | 67.6 | 238.0 | 15900 |
MATLAB中提供了两种方法可以进行聚类分析。
(1)一次聚类。它的优点是可利用clusterdata函数对样本数据进行一次聚类。它的缺点是可供用户选择的范围比较窄,不能更改距离的计算方法。
(2)分步聚类。此方法可以分以下步骤进行分步聚类。
① 找到数据集合中变量两两之间的相似性和非相似性,并用pdist函数计算变量之间的距离;
② 用linkage函数定义变量之间的连接;
③ 用cophenetic函数评价聚类信息;
④ 用cluster函数创建聚类。
下面分别介绍这两种聚类方法。
1.一次聚类
在命令窗口中输入下面的命令行:
>> X=[11978 12.5 93.5 31908
28446 30.4 202.0 23754
2501 67.2 24.8 58
1028 28.4 14.0 599
210 8.6 1.5 1147
1458 26.7 16.0 1288
635 21.1 3.6 514
32613 32.7 192.8 2385
10700 13.9 10.5 45190
92000 41.1 841.5 37370
458 35.8 8.9 168
868 27.8 11.4 405
161 17.4 2.5 129
929 36.3 11.4 640
634 26.7 11.3 447
385 34.7 2.5 200
6748 20.5 29.0 1200
2180 84.0 33.7 1200
1490 16.1 0.8 2090
4850 24.6 32.6 7450
57500 67.6 238.0 15900];
>> T = clusterdata(X,0.9)
T =
10
10
3
8
2
9
8
10
10
7
1
8
1
8
8
1
4
9
9
5
6
可见,MATLAB将数据集合分为10类。调整cutoff值,将有不同的分类。
2.分步聚类
(1)寻找变量之间的相似性
用pdist函数可以计算相似性矩阵或非相似性矩阵。在MATLAB中有多种方法可以计算距离。默认时用pdist函数计算欧氏距离,可以指定一种或多种选项(参见函数部分的内容)。在进行计算之前最好先将数据(用zscore函数)进行标准化。下面的代码返回变量之间的距离信息。
>> Y=pdist(X)
Y =
1.0e+04 *
1 ~ 13 列
1.8376 3.3230 3.3169 3.2935 3.2377 3.3380 3.6020 1.3344 8.0212 3.3766 3.3405 3.3905 3.3163
14 ~ 26 列
3.3444 3.3761 3.1150 3.2233 3.1609 2.5476 4.8255 3.5138 3.5888 3.6172 3.5116 3.6243 2.1771
27 ~ 39 列
2.7829 6.4999 3.6601 3.6135 3.6854 3.5937 3.6287 3.6637 3.1297 3.4621 3.4583 2.8681 3.0097
40 ~ 52 列
0.1570 0.2537 0.1613 0.1922 3.0202 4.5871 9.6969 0.2046 0.1670 0.2342 0.1677 0.1908 0.2121
53 ~ 65 列
0.4398 0.1186 0.2270 0.7756 5.7236 0.0985 0.0812 0.0402 3.1636 4.5628 9.8126 0.0715 0.0251
66 ~ 78 列
0.0986 0.0107 0.0422 0.0757 0.5752 0.1301 0.1561 0.7845 5.8509 0.1256 0.0763 3.2427 4.5275
79 ~ 91 列
9.8682 0.1010 0.0992 0.1019 0.0880 0.0819 0.0963 0.6538 0.1972 0.1590 0.7827 5.9160 0.1130
92 ~ 104 列
3.1175 4.4864 9.7470 0.1502 0.1062 0.1739 0.0837 0.1177 0.1528 0.5291 0.0730 0.0803 0.7034
105 ~ 117 列
5.7916 3.2033 4.5796 9.8522 0.0389 0.0257 0.0611 0.0320 0.0068 0.0402 0.6151 0.1692 0.1793
118 ~ 130 列
0.8116 5.8910 4.8088 6.8929 3.2232 3.1807 3.2531 3.1733 3.2038 3.2303 2.5893 3.0457 3.1125
131 ~ 143 列
2.8222 2.8320 8.1679 4.6172 4.5852 4.6277 4.5609 4.5861 4.6157 4.4167 4.4808 4.4073 3.8191
144 ~ 156 列
5.5210 9.8816 9.8347 9.9106 9.8202 9.8548 9.8872 9.2611 9.6833 9.7147 9.2147 4.0640 0.0474
157 ~ 169 列
0.0300 0.0667 0.0330 0.0080 0.6374 0.2008 0.2182 0.8504 5.9172 0.0759 0.0243 0.0238 0.0525
170 ~ 182 列
0.5934 0.1535 0.1796 0.8093 5.8714 0.0923 0.0570 0.0236 0.6674 0.2287 0.2369 0.8694 5.9469
183 ~ 195 列
0.0353 0.0700 0.5846 0.1372 0.1555 0.7858 5.8593 0.0351 0.6160 0.1721 0.1853 0.8174 5.8929
196 ~ 208 列
0.6441 0.2056 0.2189 0.8515 5.9234 0.4568 0.5333 0.6532 5.2838 0.1129 0.6797 5.7240 0.6326
209 ~ 210 列
5.7688 5.3324
为了便于阅读,可以用squareform函数将距离向量转换为矩阵:
>> squareform(Y)
ans =
1.0e+04 *
1 ~ 13 列
0 1.8376 3.3230 3.3169 3.2935 3.2377 3.3380 3.6020 1.3344 8.0212 3.3766 3.3405 3.3905
1.8376 0 3.5138 3.5888 3.6172 3.5116 3.6243 2.1771 2.7829 6.4999 3.6601 3.6135 3.6854
3.3230 3.5138 0 0.1570 0.2537 0.1613 0.1922 3.0202 4.5871 9.6969 0.2046 0.1670 0.2342
3.3169 3.5888 0.1570 0 0.0985 0.0812 0.0402 3.1636 4.5628 9.8126 0.0715 0.0251 0.0986
3.2935 3.6172 0.2537 0.0985 0 0.1256 0.0763 3.2427 4.5275 9.8682 0.1010 0.0992 0.1019
3.2377 3.5116 0.1613 0.0812 0.1256 0 0.1130 3.1175 4.4864 9.7470 0.1502 0.1062 0.1739
3.3380 3.6243 0.1922 0.0402 0.0763 0.1130 0 3.2033 4.5796 9.8522 0.0389 0.0257 0.0611
3.6020 2.1771 3.0202 3.1636 3.2427 3.1175 3.2033 0 4.8088 6.8929 3.2232 3.1807 3.2531
1.3344 2.7829 4.5871 4.5628 4.5275 4.4864 4.5796 4.8088 0 8.1679 4.6172 4.5852 4.6277
8.0212 6.4999 9.6969 9.8126 9.8682 9.7470 9.8522 6.8929 8.1679 0 9.8816 9.8347 9.9106
3.3766 3.6601 0.2046 0.0715 0.1010 0.1502 0.0389 3.2232 4.6172 9.8816 0 0.0474 0.0300
3.3405 3.6135 0.1670 0.0251 0.0992 0.1062 0.0257 3.1807 4.5852 9.8347 0.0474 0 0.0759
3.3905 3.6854 0.2342 0.0986 0.1019 0.1739 0.0611 3.2531 4.6277 9.9106 0.0300 0.0759 0
3.3163 3.5937 0.1677 0.0107 0.0880 0.0837 0.0320 3.1733 4.5609 9.8202 0.0667 0.0243 0.0923
3.3444 3.6287 0.1908 0.0422 0.0819 0.1177 0.0068 3.2038 4.5861 9.8548 0.0330 0.0238 0.0570
3.3761 3.6637 0.2121 0.0757 0.0963 0.1528 0.0402 3.2303 4.6157 9.8872 0.0080 0.0525 0.0236
3.1150 3.1297 0.4398 0.5752 0.6538 0.5291 0.6151 2.5893 4.4167 9.2611 0.6374 0.5934 0.6674
3.2233 3.4621 0.1186 0.1301 0.1972 0.0730 0.1692 3.0457 4.4808 9.6833 0.2008 0.1535 0.2287
3.1609 3.4583 0.2270 0.1561 0.1590 0.0803 0.1793 3.1125 4.4073 9.7147 0.2182 0.1796 0.2369
2.5476 2.8681 0.7756 0.7845 0.7827 0.7034 0.8116 2.8222 3.8191 9.2147 0.8504 0.8093 0.8694
4.8255 3.0097 5.7236 5.8509 5.9160 5.7916 5.8910 2.8320 5.5210 4.0640 5.9172 5.8714 5.9469
14 ~ 21 列
3.3163 3.3444 3.3761 3.1150 3.2233 3.1609 2.5476 4.8255
3.5937 3.6287 3.6637 3.1297 3.4621 3.4583 2.8681 3.0097
0.1677 0.1908 0.2121 0.4398 0.1186 0.2270 0.7756 5.7236
0.0107 0.0422 0.0757 0.5752 0.1301 0.1561 0.7845 5.8509
0.0880 0.0819 0.0963 0.6538 0.1972 0.1590 0.7827 5.9160
0.0837 0.1177 0.1528 0.5291 0.0730 0.0803 0.7034 5.7916
0.0320 0.0068 0.0402 0.6151 0.1692 0.1793 0.8116 5.8910
3.1733 3.2038 3.2303 2.5893 3.0457 3.1125 2.8222 2.8320
4.5609 4.5861 4.6157 4.4167 4.4808 4.4073 3.8191 5.5210
9.8202 9.8548 9.8872 9.2611 9.6833 9.7147 9.2147 4.0640
0.0667 0.0330 0.0080 0.6374 0.2008 0.2182 0.8504 5.9172
0.0243 0.0238 0.0525 0.5934 0.1535 0.1796 0.8093 5.8714
0.0923 0.0570 0.0236 0.6674 0.2287 0.2369 0.8694 5.9469
0 0.0353 0.0700 0.5846 0.1372 0.1555 0.7858 5.8593
0.0353 0 0.0351 0.6160 0.1721 0.1853 0.8174 5.8929
0.0700 0.0351 0 0.6441 0.2056 0.2189 0.8515 5.9234
0.5846 0.6160 0.6441 0 0.4568 0.5333 0.6532 5.2838
0.1372 0.1721 0.2056 0.4568 0 0.1129 0.6797 5.7240
0.1555 0.1853 0.2189 0.5333 0.1129 0 0.6326 5.7688
0.7858 0.8174 0.8515 0.6532 0.6797 0.6326 0 5.3324
5.8593 5.8929 5.9234 5.2838 5.7240 5.7688 5.3324 0
(2)定义变量之间的连接
将数据中变量之间的近似性计算出来以后,使用linkage函数就可以决定数据集合中哪些变量可以被归入某类中。该函数利用pdist函数生成的距离信息,连接相近的成对变量,形成二分类。然后,linkage函数将新生成的类与其他变量相连接,生成更大的类,直到原始数据中所有的变量被连接到系统聚类树中。
>> Z=linkage(Y)
Z =
1.0e+04 *
0.0007 0.0015 0.0068
0.0011 0.0016 0.0080
0.0004 0.0014 0.0107
0.0013 0.0023 0.0236
0.0012 0.0022 0.0238
0.0024 0.0026 0.0243
0.0025 0.0027 0.0330
0.0006 0.0018 0.0730
0.0005 0.0028 0.0763
0.0019 0.0029 0.0803
0.0030 0.0031 0.0812
0.0003 0.0032 0.1186
0.0017 0.0033 0.4398
0.0020 0.0034 0.6326
0.0001 0.0009 1.3344
0.0002 0.0036 1.8376
0.0008 0.0037 2.1771
0.0035 0.0038 2.5476
0.0021 0.0039 2.8320
0.0010 0.0040 4.0640
在输出信息中,每一行代表一次连接,前两列代表连接的变量,第3列代表连接对象之间的距离。例如,首先将国家7和国家15相连接,形成新类22,它们之间的距离是68;然后将国家16和国家11相连,形成新类23,它们之间的距离为80,比前一组略大;接着,将国家4和国家14相连,形成新类24;再将新类23与国家13相连,形成新类25;依此类推,最后将新类40和国家10相连,形成整个系统。整个输出信息展示了聚类的过程。
(3)评价聚类信息
连接变量生成聚类树以后,用户可能会想修改聚类树或希望了解变量连接的更多信息。下面分别介绍具体实现方法。
① 修改聚类树:衡量通过linkege函数生成的聚类信息的有效性的一种方法是将它与pdist函数生成的原始相似数据对比。如果有效,则二者强相关。可以使用cophenet函数比较这两组值,并且计算它们的相关性,返回cophenetic相关系数。该值越接近于1,表示聚类效果越好。
可以通过cophenetic相关系数来比较用不同距离计算方法或聚类算法对同一套数据得到的聚类结果。对于本例:
>> C=cophenet(Z,Y)
C =
0.9393
将pdist函数中距离的计算方法分别指定为“Mahal”“SEuclid”和“CityBlock”,重新使用pdist函数计算以后,用cophenet函数计算,得到的C值分别等于0.7486、0.8018和0.9291。均小于距离为欧氏距离时的计算结果,所以使用默认设置时的效果最佳。
② 了解与聚类连接相关的更多信息:决定数据集合中聚类的方法之一是比较聚类树中每一个连接的长度与相邻次一级连接的长度。如果二者相近,则表示在此水平上变量之间是相似的,这些连接被认为具有较高水平的连续性;如果二者不相近,则表示聚类树中变量之间不相似,则称该连接与周围是不连续的。cluster函数用不连续性的度量决定在哪里对数据进行分类。
下面的命令行生成聚类树:
>> dendrogram(Z)
ans =
151.0055
153.0055
154.0055
155.0055
156.0055
157.0055
158.0055
159.0055
160.0055
161.0055
162.0055
163.0055
164.0055
165.0055
166.0055
167.0055
168.0055
169.0055
170.0055
171.0055
图5-1中的纵轴为变量之间的距离,横轴为聚类变量。系统聚类树从左到右反映了聚类的先后次序。该图可以与linkage函数的计算结果相对应。
图5-1 本问题的冰柱图
在系统聚类树中,每一个连接的相对连续性可用不连续性系数定量表达,该函数用于比较某个连接的长度与相邻连接长度的均值。若该变量与周围变量连续,则不连续性系数较低,反之则反。
用inconsistent函数可以生成聚类树中每一个连接的不连续系数表。该函数可以比较每一个连接与以下两个水平的连接,这里称为比较的深度。用该函数可以指定其他深度。聚类树底部的对象被称为叶节点,不连续系数为0。
>> I=inconsistent(Z)
I =
1.0e+004*
0.0068 0 0.0001 0
0.0080 0 0.0001 0
0.0107 0 0.0001 0
0.0158 0.0110 0.0002 0.0001
0.0153 0.0120 0.0002 0.0001
0.0196 0.0077 0.0003 0.0001
0.0270 0.0053 0.0003 0.0001
0.0730 0 0.0001 0
0.0546 0.0306 0.0002 0.0001
0.0766 0.0052 0.0002 0.0001
0.0793 0.0026 0.0003 0.0001
0.0999 0.0265 0.0002 0.0001
0.2792 0.2271 0.0002 0.0001
0.5362 0.1363 0.0002 0.0001
1.3344 0 0.0001 0
1.5860 0.3559 0.0002 0.0001
2.0074 0.2401 0.0002 0.0001
1.7858 1.0157 0.0003 0.0001
2.6898 0.2011 0.0002 0.0001
3.4480 0.8711 0.0002 0.0001
在矩阵中,第1列为所有连接长度的均值,第2列为所有连接长度的标准差,第3列为计算中所包含的连接数,第4列为不连续系数。该输出信息可以与linkage函数的输出信息进行对比,第1行代表国家7和国家15相连,因为二者均为叶节点,所以不连续系数为0。
(4)创建聚类
创建二分聚类树以后,可以用cluster函数将系统分成更大的类。有以下两种方式可以创建聚类:
① 找到原始数据集合的自然分界:如果用cluster函数对数据进行分类,指定一个不连续系数0.9作为cutoff参数的值,则cluster函数将所有变量分为10类,如下所示。
>> T=cluster(Z, 0.9)
T =
10
10
3
8
2
9
8
10
10
7
1
8
1
8
8
1
4
9
9
5
6
如果选择不连续系数0.7作为cutoff参数的值,则所有数据分为16类。
>> T=cluster(Z, 0.7)
T =
16
9
6
3
4
15
14
10
16
12
13
2
1
3
14
13
7
15
5
8
11
② 通过指定类数进行聚类:将所有变量归为3类,将cutoff参数的值设置为3,用cluster函数进行分类。
>> T=cluster(Z,3)
T =
2
2
2
2
2
2
2
2
2
3
2
2
2
2
2
2
2
2
2
2
1