前言
1992年6月，联合国在里约热内卢召开的“环境与发展大会”，通过了以可持续发展为核心的《里约环境与发展宣言》、《21世纪议程》等文件。随后，中国政府编制了《中国21世纪人口、资源、环境与发展白皮书》，首次把可持续发展战略纳入我国经济和社会发展的长远规划。1997年的中共十五大把可持续发展战略确定为我国“现代化建设中必须实施”的战略。2002年中共十六大把“可持续发展能力不断增强”作为全面建设小康社会的目标之一。
而即使是在已经达成小康社会的今天，可持续发展依然是中国需要遵守的基本科学发展观。于是我们提出问题：能否运用R语言的多元统计分析，对当今某一省的可持续发展程度做出客观评价。
本文选择广东省做分析，介绍了关于可持续发展的指标选取理由。从现有文献来看，运用R软件进行实证分析的较少，多数还是用文字描述后提出建议。因此本文将根据《广东省统计年鉴》的统计数据，采用对数据实证分析的方法，构建指标体系，通过因子分析和聚类分析等方法评价对广东的可持续发展水平进行分析，最后给出一些提议，用来弥补现有文献的不足和缺陷。
1 指标选取
从可持续发展的角度看，我们选取的指标一定要详细、全面的反映出可持续发展的现状、水平、趋势。因此我们寻找一组具有非常典型代表意义的指标，这些指标将符合一些标准。
第一是相关性。作为反映某个省的可持续发展指标，各指标之间可以长期保持有序、协调的关系，以便于人们可以清楚的掌握到该省的持续发展程度。这就要求指标体系里每个指标间有着紧密的联系。
第二是可操作性。现如今，环保指标，经济指标等等，为了全面的反映现实状况，一般数量都非常惊人，少也有上百个。而对于可持续发展的指标来说，数量过多、过于复杂的一些指标往往不好操作，而且有些指标难以获取精确的数据。因此，我们在建立指标体系的时候，要以简便为前提，选择一些容易计算，方便获取精确数据且能很好的反映该省水平的指标，从而使该分析具有可操作性。
第三是动态性。可持续发展本身就是长期持续的过程，不同的时期指标内涵也一定会有所不同。这就需要我们选取的指标不能只反映某一个阶段的问题，还要有一定的弹性，能够反映出不同阶段该省的发展特征。
基于上面三条原则，本文选取了地区生产总值(单位：亿元)、人均地区生产总值(单位：万元)、固定资产投资额(单位：亿元)、进出口总额(单位：亿美元)、第三产业 比重 (%)、年末城镇失业率 (%)、中等职业教育在校学生数 (单位：万人)、医疗保险参保人数 (单位：万人)、医院床位数(单位：张)、职业医师数(单位：人)、城镇污水处理率(%)、工业废水排放总量 (单位：亿吨)、工业烟尘排放总量 (单位：万吨)以上十三个指标。
从因子分析角度来看，数据中多个变量之间要有较强的线性相关关系。这个将会在下面的分析中提到。
2 基于R软件的实证分析
2.1数据来源
本文的所有数据均来自广东统计局年鉴，选取了上述十三个指标反映前一年的可持续发展状况，本文用x1-x13表示这十三个指标。数据详情见附录一。
2.2分析方法
主成分分析的基本思想：
主成分分析是设法将原来众多具有一定相关性（比如P个指标），重新组合成一组新的互相无关的综合指标来代替原来的指标。主成分分析，是考察多个变量间相关性一种多元统计方法，研究如何通过少数几个主成分来揭示多个变量间的内部结构，即从原始变量中导出少数几个主成分，使它们尽可能多地保留原始变量的信息，且彼此间互不相关.通常数学上的处理就是将原来P个指标作线性组合，作为新的综合指标。
因子分析的基本思想：
因子分析是通过降维的手段将复杂的许多变量归为少数的几个综合因子的一种分析方法，也就是说在具有一定相关性的原变量中提取公因子，通过构建因子模型用来描述原变量之间的相关关系，目前因子分析在心理学、社会学、经济学等学科中都有重要的应用。
2.3 因子分析使用性检验
我们采用KMO和Bartlett检验两种方法来检验各变量间有没有相关关系。利用R软件对选取的十三个指标进行因子分析，所得到的的总KMO值和各个变量的KMO值见图1。
图一 KMO检验
由图2可得到，总体KMO检验值是0.66，超过了0.6这个临界值，比较适合做因子分析。
R软件得到的Bartlett球形检验结果见图2
图2 Bartlett球形检验
Bartlett球形检验的结果说明p值远小于0.05，适合做因子分析。
2.4 主成分分析
2.4.1 计算相关矩阵
相关矩阵见图3
图3 相关矩阵
2.4.2 求相关矩阵的特征值和主成分负荷
特征值见图4
图4 特征值
其中特征值是开根号后的结果。standard deviation，proportion of variance，cumulative proportion分别代表标准误，方差比率，累计贡献率。
主成分载荷见图5
图5 主成分载荷表
2.4.3 确定主成分
按照累计方差贡献率大于80%的原则，选定了四个主成分，其累计方差贡献率为85.8%，所以取m=4。从碎石图也可以看出m取4比较合适。碎石图见图6。
图6 碎石图
2.4.4 主成分得分及结果分析
主成分得分见图7
图7 主成分得分
结果分析：由主成分载荷矩阵可以看出，主成分Comp.1在x1（地区生产总值）,x2（人均地区生产总值）,x3（固定资产投资额）,x5（第三产业比重）,x8（医疗保险参保人数）,x9（医院床位数）,x10（职业医生数）,x12（工业废水排放量）上的载荷值都很大，可视为反映出了经济和医疗方面的可持续发展水平；Comp.2在x7（中等职业教育在校学生数）上的载荷值很大，可视为反映出了教育的可持续发展水平；Comp.3在x4（进出口总额），x13（工业烟尘排放总量）上载荷值很大，可视为反映出了外贸和工业的可持续发展水平；Comp.4在x6（年末城镇失业率），x11（城镇污水处理率）上的载荷值很大，可视为反映出了就业和环保方面的可持续发展水平。
这样，有了各个主成分的解释，再结合各个地区在四个主成分上的得分和综合得分，就可以对各个地区进行评价了。最后，由加权法估计综合得分，以各主成分的方差贡献率占四个主成分总方差贡献率的比重作为权重进行加权汇总，得到各个地区的综合得分。
各个地区的综合得分见图8。（PC为综合得分）
图8 综合得分
从排名可以看出来作为中国的一线城市，广州和深圳的综合得分远远超过了其他城市，深圳虽然在经济和医疗方面不如广州，但在就业、教育、环保方面会比广州强。而其他城市也有些有自己的可持续发展优点，比如梅州在外贸和工业方面比广州都强。云浮和河源的各方面指标评分都比较低，可持续发展能力较弱，地方政府可以着重照顾一下，带动这些地区持续发展。
2.5 因子分析
2.5.1 选择因子个数
由上面的主成分分析法我们得知了前四个因子就能概括大部分信息，这次我们用主成分法因子分析再看看累计方差贡献率，见图9。
图9 旋转前后因子方差贡献
可见旋转前后依然是前四个因子的累计方差贡献率超过80%，所以取前四个因子就可以很好的覆盖原指标。
2.5.2 因子载荷矩阵
由R给出的旋转前后因子载荷矩阵（图10）可以得出，Factor1在x1（地区生产总值），x2（人均地区生产总值），x3（固定资产投资额），x4（进出口总额），x5（第三产业比重）这四个指标上的载荷很大，可视为反映了可持续发展的经济状况；Factor2在x6（年末城镇失业率）上的载荷很大，可视为反映了可持续发展的就业情况；Factor3在x7（中等职业教育在校学生），x8（医疗保险参保人数），x9（医院床位数），x10（职业医生数）上的载荷很大，可视为反映了可持续发展的教育和医疗情况；Factor4在x11（城镇污水处理率），x12（工业废水排放总量），x13（工业烟尘排放总量）的载荷很大，可视为反映了可持续发展的环保情况。
图10 旋转前后因子载荷
2.5.3 因子排名分析
因子得分及排名见图11和图12，由此我们可以看出，深圳在可持续发展经济方面，远远超过其他城市，稳居第一，但其他指标发展的却并不统一，尤其是环境和医疗教育方面，深圳虽然综合评分的排名位列第二，但也是全靠经济，属于非常典型的经济发展过快而其他方面跟不上。相比之下广州的各指标发展都很均衡，只在医疗教育方面有些弱，所以才位列第一。所以想要全面的可持续发展就必须保证各个指标全面发展，不能顾此失彼。位列第三的是佛山，作为非一线城市，虽然评分不如广州深圳，但佛山也有其突出优点，即医疗教育远高于其他城市，但经济发展不太景气，所以地方政府可以加大对佛山的投资，让其经济发展也一并跟上去。湛江的环保评分是最突出的，事实上，湛江也确实被称为广东的花园城市，因此在其他指标远不如广州深圳的情况下，依然位列第四，由此可见环保对可持续发展的重要性，地方政府可以向湛江输送人才来建设家乡。东莞在各指标的综合发展也很均衡，属于是效仿广州而成功的一个非常好的例子，所以位列第五。接下里的茂名，在地理位置上紧挨着湛江，所以茂名的环保评分也很高。
综上所述，可持续发展一定要全面，均衡的发展，不能只看重某一个指标，尤其是不能牺牲环境来让经济发展，走西方的老路。
图11 旋转后因子得分
图12 旋转后因子得分排序
我们以因子1和因子2作为坐标轴，绘制各城市的得分图和信息重叠图，见图13和图14。由此可以看出广州深圳这两个一线城市卓越的经济实力，各个变量在广州和深圳这两个城市在信息重叠图里的反映也很强烈。
图13 旋转后因子得分图
图14 信息重叠图
2.6 聚类分析
在因子分析结束后，我们对各变量进行聚类分析，对广东省各城市的可持续发展情况由高到低进行分类和排序，目的是将分析的结果更加一目了然。系统聚类有最短距离法、最长距离法、中间距离法、类平均法、重心法、ward法六种，这里我们选取类平均法，可以综合利用所有变量的信息。见图15。
图15 聚类分析图
由此我们可以将所有城市分为四类，
|
类别 |
城市名称 |
|
第1类 |
广州 |
|
第2类 |
深圳 |
|
第3类 |
珠海，佛山，中山，东莞， |
|
第4类 |
汕头，韶关，河源，梅州，惠州，汕尾，江门，阳江，湛江，茂名，肇庆，清远，潮州，揭阳，云浮 |
我们可以看到，等级划分和现实情况基本一致。广州作为全面发展的一线城市位居第一类；深圳没有照顾好经济外的指标位居第二类；第三类属于有自己的突出优点或发展的比较全面的城市；第四类属于优点突出不明显、各方面发展都中等偏下。
3 改善建议
3.1 重新对城市进行规划，使之布局合理
我们可以看到在广东的所有城市中，只有湛江、茂名这几个少数城市在环保指标上比较优秀，所以可以通过政府出台政策，控制城市土地价格，约束城市恶性扩张土地。在可持续发展的过程中，不仅要考虑到经济的发展，还要考虑到环境的保护。比如重工业区在扩张的时候，一定要在下风口处布局，注重河流、泥沙的流向，考虑到河流的自净能力。
3.2 探索“绿色经济”
加快北部生态发展区绿色发展步伐，完善生态补偿机制，支持各市壮大县域、镇域经济，因地制宜发展水经济、清洁能源、绿色矿业等，积极发展农产品种养和深加工、乡村旅游经济，打造生态经济发展新标杆。广东省深化新一轮省内对口帮扶，引导珠三角产业向粤东粤西粤北地区梯度转移，推动各功能区深化产业统筹、项目对接。粤东粤西粤北各市要下大决心优化营商环境，大抓招商引资，增强产业承接能力，打造珠三角产业拓展首选地和先进生产力延伸区。
3.3 使地方优势更加突出
在云浮市的中药材种子种苗繁育基地，将与中国中医科学院合作共建世界一流的药用资源种质库，打造区域性的种子种苗繁育中心。截至目前，云浮市建设种子种苗的繁育基地总共有15个，种子种苗的繁育量达到1亿余株，有望建成世界一流的药用资源种质库。此外，广东北部的梅州市、韶关市、清远市等地区的“绿色经济”也取得了骄人的成绩。在这其中，清远市围绕清远鸡、英德红茶、连州菜心、丝苗米、麻竹笋等特色农产品，加速形成五个百亿产值产业。除了传统产业升级外，清远、韶关、河源、梅州以及云浮五市还瞄准发展新兴产业。
上述城市都是在分析中不及广州深圳的城市，但它们也可以通过扩大自己的优势得到属于自己的可持续发展，以“一核一带一区”区域发展格局，实现更高水平、更高质量的区域协调发展，广东全省努力把短板变成“潜力板”，加快形成主体功能明显、优势互补的区域经济布局，奋力在全面建设社会主义现代化国家新征程中走在全国前列。
4 结语
本文通过分析广东省统计年鉴的数据，在阐述了选取可持续发展指标的原则后，建立指标体系，利用主成分分析法，因子分析法，聚类分析法，对广东省这一年的可持续发展情况作出客观评价。最后提出了一些改善的建议。因为笔者能力有限，难免会有疏忽和失误，希望后人能对广东省可持续发展做出进一步分析，也希望本文的建议能对解决或改善广东省的可持续发展问题有所帮助。
参考文献
- 焦子涵. 基于R语言因子分析和聚类分析的市政建设水平综合评价研究[J]. 福建建材, 2022(10):5.
- 祖来克孜米吉提. 基于R语言因子分析法的区域经济发展水平综合评价[J]. 当代经济, 2021.
- 林双双. 中国房地产价格主要影响因素研究——基于R语言因子分析[J]. 福建质量管理, 2019.
- 张翠娟, 冯学军, 盛敏. 因子分析开发步骤及R语言程序代码实现[J]. 安庆师范学院学报：自然科学版, 2013, 19(2):4.
- 方匡南. 基于数据挖掘的分类和聚类算法研究及R语言实现[D]. 暨南大学, 2008.
附录一 实验数据
|
x1 |
x2 |
x3 |
x4 |
x5 |
x6 |
x7 | |
|
广州 |
13551.21 |
105908.94 |
3758.39 |
1171.67 |
63.59 |
2.41 |
24.01 |
|
深圳 |
12950.06 |
123247.05 |
2314.43 |
4668.03 |
55.65 |
2.42 |
3.21 |
|
珠海 |
1503.76 |
95471.03 |
787.62 |
456.8 |
45.78 |
2.34 |
2.27 |
|
汕头 |
1425.01 |
26230.79 |
611.92 |
88.02 |
42.41 |
2.43 |
10.97 |
|
佛山 |
6613.02 |
91259.42 |
2128.33 |
610.58 |
35.83 |
2.42 |
7.98 |
|
韶关 |
906.48 |
31702.17 |
548.48 |
20.44 |
44.54 |
2.75 |
4.02 |
|
河源 |
615.26 |
20536.38 |
278.59 |
29.26 |
38.57 |
2.5 |
3.89 |
|
梅州 |
744.75 |
17396.3 |
230.14 |
15.03 |
42.76 |
2.4 |
7.07 |
|
惠州 |
2367.55 |
50873.48 |
1208.68 |
494.94 |
36.5 |
2.35 |
6.99 |
|
汕尾 |
610.41 |
20607.96 |
391.56 |
28.42 |
37.18 |
2.62 |
5.7 |
|
东莞 |
5010.17 |
60556.5 |
1180.35 |
1445.17 |
52.21 |
2.3 |
4.73 |
|
中山 |
2441.04 |
77526.66 |
893.43 |
335.23 |
42 |
2.3 |
2.55 |
|
江门 |
1880.39 |
42028.43 |
850.41 |
187.72 |
40.95 |
2.35 |
5.23 |
|
阳江 |
887.03 |
36095.67 |
483.67 |
22.23 |
34.06 |
2.5 |
3.04 |
|
湛江 |
1860.22 |
26240.21 |
572.28 |
47 |
40.43 |
2.75 |
12.82 |
|
茂名 |
1936.18 |
32677.57 |
427.37 |
10.41 |
41.45 |
2.75 |
11.91 |
|
肇庆 |
1462.35 |
36864.27 |
852.6 |
63.52 |
37.91 |
2.42 |
6.91 |
|
清远 |
1025.03 |
27319.62 |
437.95 |
45.37 |
44.46 |
2.75 |
4.98 |
|
潮州 |
706.65 |
26251.62 |
224.16 |
42.31 |
37.99 |
2.35 |
2.38 |
|
揭阳 |
1396.79 |
23532.3 |
663.51 |
42.75 |
28.44 |
2.42 |
14.51 |
|
云浮 |
530.29 |
22115.33 |
463.66 |
14.57 |
34.72 |
2.75 |
4.42 |
|
x8 |
x9 |
x10 |
x11 |
x12 |
x13 | ||
|
广州 |
768.22 |
70649 |
37101 |
82.73 |
2.27 |
1.42 | |
|
深圳 |
1138.74 |
28065 |
23973 |
96.1 |
1.38 |
0.52 | |
|
珠海 |
150.09 |
7439 |
4733 |
86.55 |
0.55 |
1.12 | |
|
汕头 |
491.16 |
13527 |
7894 |
90.23 |
0.51 |
0.49 | |
|
佛山 |
448.19 |
25686 |
13171 |
99.04 |
2.04 |
4.24 | |
|
韶关 |
87.45 |
13510 |
6058 |
81.88 |
1.02 |
0.66 | |
|
河源 |
63.95 |
9228 |
4042 |
89.39 |
0.19 |
0.52 | |
|
梅州 |
459.13 |
12660 |
7399 |
53.34 |
0.57 |
0.57 | |
|
惠州 |
392.37 |
17231 |
8564 |
92.06 |
0.83 |
2.13 | |
|
汕尾 |
289.6 |
6766 |
4345 |
85.57 |
0.19 |
1.09 | |
|
东莞 |
616.86 |
24617 |
12849 |
95.11 |
2.69 |
1.51 | |
|
中山 |
250.68 |
11336 |
5635 |
90.7 |
0.86 |
1.24 | |
|
江门 |
379.85 |
15627 |
7409 |
88.15 |
1.44 |
1.02 | |
|
阳江 |
247.37 |
8442 |
3381 |
83.5 |
0.2 |
1.17 | |
|
湛江 |
693.09 |
24212 |
9587 |
96.64 |
0.82 |
1.01 | |
|
茂名 |
71.87 |
19748 |
9624 |
85.81 |
0.52 |
1.04 | |
|
肇庆 |
394.63 |
11869 |
5176 |
85.39 |
1.26 |
2.48 | |
|
清远 |
400.01 |
11823 |
6381 |
83.57 |
0.41 |
3.13 | |
|
潮州 |
251.34 |
5602 |
3636 |
86.2 |
0.3 |
0.33 | |
|
揭阳 |
566.22 |
10812 |
7941 |
80.82 |
0.37 |
0.26 | |
|
云浮 |
261.01 |
6425 |
3208 |
98.56 |
0.17 |
0.81 |
附录二 R语言代码
date=read.table("clipboard",header=T) #从剪贴板里读取数据
library(psych) #KMO检验的程序包
library(mvstats) #主成分分析和因子分析要用的程序包
KMO(date) #KMO检验
bartlett.test(date) #bartlett检验
cor(date) #相关系数矩阵
PCA=princomp(date,cor=T) #主成分分析
summary(PCA) #显示标准误，方差贡献率等因素
PCA$loadings #主成分载荷
screeplot(PCA,type="lines") #碎石图
PCA$scores[,1:4] #主成分得分
princomp.rank(PCA,m=4) #主成分排名
princomp.rank(PCA,m=4,plot=T) #主成分作图
FA0=factpc(date,4) #因子分析（不旋转）
FA0 #显示方差和贡献率（不旋转）
Fa1=factpc(date,4,rot="varimax") #因子旋转
Fa1$Vars #旋转后的方差
Fa1$loadings #旋转后的因子载荷
Fa1$scores #旋转后的因子得分
Fa1$Rank #旋转后的因子排名
plot.text(Fa1$scores) #因子得分图
biplot(Fa1$scores,Fa1$loading) #信息重叠图
hc=H.clust(date,"euclidean","average",plot=T) #类平均法的聚类分析
cutree(hc,4) #分为四类
Comments NOTHING