【毕业学位论文】（Word原稿）北京市水资源短缺风险综合评价——基于因子分析优化模型-统计教育学

第 1 页 /共 27 页 北京市水资源 短缺风险 综合评价 基于因子分析 优化模型 华南师范大学 摘要 近年来 ，我国， 特别是北方地区 的 水资源短缺问题日趋严重，水资源 短缺 成为 了 焦点话题。 作为我国 的 首都，北京市 也 是一座严重缺水的特大城市，已连续干旱了十多年。 为了更 好的了解我国 首都 北京的水资源的基本情况 以及找出其水资源短缺的主要原因 ， 对其作 出 综合评价 以 及给当地政府一些相关的 参考 建议 。 本论文从北京市 2010 年度 统计年鉴找到 2001 年到 2009 年北京市影响水资源的 相关指标的数据， 主要应用 因子分析法对其 进行分析 ， 找出影响北京市水资源短缺的主要风险因子 以及综合风险评价因子 ，同时 在此基础上 建立 了因子分析优化 预测模型， 对 9 年来水资源的情况进行判别； 并 用一次指数平滑法对该市 2010 年的有关数据进行 估计 ， 以此预测 2010 年水资源的好坏 ，并进一步对多年来影响北京市水资源问题的原因进行分析 以及提出解决的建议 。 关键词： 指标体系， 因子分析 ， 一次指数平滑 估计，因子分析优化 模型。 第 2 页 /共 27 页 图片来源： 2009年北京市水资源公报 正文 一、 问题 阐述 北京市是我国一座严重缺水的特大城市， 已连续干旱了十多年 ；此外，降雨量又没有增加反而在逐年下降 (如图 1以致 其人均占有量仅为全国平均水平的 1/8，为世界人均水平的 1/30，属于重度缺水地区，水资源短缺的问题已成为制约其发展的主要因素。尽管政府采取了如南水北调工程建设、产业结构调整等一系列的措施，但气候变化和经济社会 在 不断发展，水资源短缺风险始终存在，且日趋严峻。 北 京市水资源 公报中非常详细地记录了北京市水资源各年的状况，包括其概述、水资源、水资源利用、水质和重要水事。要找到水资源短缺的主要风险因子，首先要知道影响水资源的因素有哪些，再从分析众多因素对水资源影响的大小程度来判定该因素是否是主要风险 因子。 近几年来北京市水资源的短缺风险程度究竟是怎样的 ？ 又该 怎样 作出综合的 评价？ 国家统计局对 2010 年水资源的有关研究工作结果还没有出来，能否 借助前人的经验建立主成分分析模型 1或者因子分析 模型 对北京市水资源短缺风险进行年度性 研究，预测 并反映出 2010 年的用水 风险程度 呢？ 050100150200降水量（ 2 3 4 5 6 7 8 9 10 11 12 月2009年 2008年 多年平均图 108、 09 年降雨量与多年平均对比 第 3 页 /共 27 页 二、 数据描述 量 选取 在考虑水资源短缺指标时，要考虑水资源质量和数量的指标，也要考虑水资源的开发利用程度和水资源开发利用效率方面的指标。 综合分析，本论文考虑了 影响水 资源 短缺风险 的水资源压力指标、水资源消耗指标、污水响应指标 和 经济 效益指标 ，这四类指标涵盖了供水、用水、气象、经济等方面对水资源的影响 。而水资源压力 指标 方面的变量 有降雨量、地表水资源量、地下水资源量、水资源总量、总供水量及 人均水资源量；水资源消耗指标方面的变量 有工业用水量、 人均年 生活用水量、农业用水量和环境用水量；污水响应指标 方面的变量 有污水处理率、污水处理能力、工业废水排放量和生活污水排放量；经济效益指标 方面的变量 有万元 耗量和万元地区生产总值水耗下降率。 (原始数据看附表 I) 标体系建立 通过 所选取的变量，建立能解决问题的指标体系至关重要 。 在运用因子分析法 对变量进行分析时，由于 所选 的 变量 是 要求 正向化的， 所以为了方便 后面的 正向化 ，我们 对以上 选取 的非正向化的 变量进行如下处理： 工业用水率 =工业用水量 /总供水量， 农业用水率 =农业用水量 /总供水量， 环境用水率 =环境用水量 /总供水量， 人均年生活用水率 =人均年生活用水量 /人均水资源量， 工业废水排放率 =工业废水排放量 /(工业废水排放量 +生活污水排放量 )， 生活污水排放率 =生活污水排放量 /(工业废水排放量 +生活污水排放量 )。 (处理后的数据看附表 经过上面的数据处理，我们根据所选变量 建立以下指标体系 (正向化的指标处理，请参阅 第 4 页 /共 27 页 表 2量明细表 选取目标 选取准则 选取变量名称 正向化分析 变量名称 北京市区域水资源短缺风险 水资源压力指标 降雨量 (单位 :毫米 ) 地表水资源量 (单位 :亿立方米 ) 地下水资源量 (单位 :亿立方米 ) 水资源总量 (单位 :亿立方米 ) 总供水量 (单位 :亿立方米 ) 人均水资源量 (单位 :立方米 ) 降雨量 ( 1X ) 地表水资源量 ( 2X ) 地下水资源量 ( 3X ) 水资源总量 ( 4X ) 总供水量 ( 5X ) 人均水资源量 ( 6X ) 水资源消耗指标 工业用水率 (%) 人均年生活用水率 (%) 农业用水率 (%) 环境用水率 (%) 工业节水比 ( 7X ) 人均年生活节水比 ( 8X ) 农业节水比 ( 9X ) 环境节水 比 ( 10X ) 污水响应指标 污水处理率 (%) 污水处理能力 (单位 :万立方米 /日 ) 工业废水排放 率 (%) 生活污水排放 率 (%) 污水处理率 ( 1X ) 污水处理能力 ( 12X ) 工业废水 节排比 ( 13X ) 生活污水 节排比 ( 14X ) 经济效益指标 万元 耗量 (单位 :立方米 ) 万元地区生产总值水耗下降率 (%) 单位水产 益 ( 15X ) 万元地区生产总值水耗下降率 ( 16X ) 型假设 为了加强我们所建立的因子分析优化模型的严密性，我们提出以下假定： (1)假设近十年北京市的自然条件是相对平稳的，即没有重大的水资源破坏 ， 没有不规律的 降雨、降雪 等； (2)假设 所选择的分析变量在短时间内都是相对稳定的且一致有效的，即观测值存在波动， 但 不会出现极端值 ，同时也不是估计的 ，从而排除了数据突变 和无效估计 带来的影响； (3)假设模型中 所有的因子都满足正态分布 ，即 2, 其中16,2,1 i ； (4)假设对 2010年各个指标的估计所产生的误差是可以忽略不计的，即 第 5 页 /共 27 页 估计的误差项是收敛于零的。 三、 模型建立 子分析优化模型理论 因子分析优化模型：有 p 维的可观测随机向量 ),(21 ，均值为 ),( 21 ，协差阵为 )( ,对正整数 ，求：称之为因子载荷阵的 ，称之为因子的随机向量 ),(21 ，称之为误差的随机向量 ),(21 p ，使 (),1,:( 2211 即 0),( ， 0)( )( ， 0)( E ， ()( 达到最大 (方阵的迹 )， (B 的列数 较靠近简单结构 )/(|m ， (式 (的关系构成正交因子优化模型。 用式 (易证明：)/( 2/1是变量 X 与因子 f 的相关阵。 注 1: 简单结构：)/( 2/1每行中有一元素绝对值靠近 1，其余元素绝对值小额到中等。 pi ,2,11 2 (B 第 j 列元素的平方和 )称为因子方差贡献，它是衡量因子 信息程度的指标； )(1 称为因子累计方差贡献，它是衡量因子 f 解释变量 X 信息程度的指标； )(称为因子 f 解释变量)(1 称为因子 f 解释变量 X 的误差信息。 引理 1 设 的特征值为 1 、p， 1 p 0， ),()(1 ， （ 这里 有 ,),(1 PP g p ） , ),( 2/112/110 kk ， 第 6 页 /共 27 页 )(,),( 2/112/110 (前 k 个标准化主成分 )， 0 ),( 1 pk )(,),( 1 ，则 000 、 足式 (式 (且 100 )( 。 定理 1 如果 0B 满足条件 (则 0 ， )( 00 ， 0 是因子分析优化模型的精确解，且1)(m a x ( , 000同引理 1)， 0f 称为初始因子， 0B 称为初始因子载荷阵。 当变量与初始因子的相关阵)/( 2/10 没有达到简单结构时，初始因子 0B 难以和实际问题对应，此时可以通过某个正交阵 ，使0f 和 0B )( 有鲜明的实际意义，即变量与旋转后因子的相关阵 )/( 2/1较靠近简单结构。这里有： 定理 2 设 是使 0B 方差最大化的正交旋转矩阵，如果 0B 满足条件 (则 , , 0 是因子分析优化模型的精确解， 且1)(m a x ( , 000同引理 1)， f 称为旋转后因子 ( 0f 的旋转 )， B 称为旋转后因子载荷阵。 标的正向化 考虑到 指标 工业用水率 、 人均年生活用水率 、 农业用水率 、 环境用水率 、工业废水排放率 、 生活污水排放率 、 水耗万元 是逆向 或者反向指标，而并 不是正向 指标 。 因此 ， 要对 这几个指标分别 进行正向化 ： 分别对 工业用水率 、人均年生活用水率 、 农业用水率 、 环境用水率 取倒数 ，并将其命名为 7X 工业 第 7 页 /共 27 页 节水比 、 8 9 10 13 14X 生活污水节排比 ，对 水耗万元 倒数后再乘以10000， 并将其命名为 效益单位水产 G D X。 （正向化后的数据看附表 关阵 的特征值 ， 21(其中 ) 表 3相关阵特征值矩阵 （ of 序号 特征值 前后特征值差 方差贡献率 累计贡献率 1 0 1 2 3 4 5 6 相关阵特征值计算编程看附录 从相关系数矩阵的特征值方差贡献率可以看到，当 时，此时的累计方差贡献率已经达到 而当 时，累计方差贡献率已达 因此，若 选 用初始因子构建因子分析优化模型，则可参考选用初始因子个数为2m 或者 3m 。 此外，我们也注意到 当 9m 时， 0d (其中 16,10,9 d )，我们查阅了相关的资料，资料显示由于所建立的指标体系中的变量数 p 大于样品量 n ，以致相关系数矩阵的协差阵在 之后就都是零， 这 正 是小样本 引起 的结果，说明后续的模型建立需要考虑小样本的影响。 第 8 页 /共 27 页 子载荷阵 表 3初始因子载荷阵 0 变量 因子 1 因子 2 因子 3 因子 4 因子 5 初始因子载荷阵计算编程看附录 从显著性 (，当因子取到 3m 时，往后就 再也 没有显著的因子了，说明选用初始因子进行因子分析最多只能选定三个因子 对所有变量进行分析，这在一定程度上会造成因子的损失以致变量的可解析度下降，同时因子分析模型的优化程度也随着下降 。尽管达到了降维的目的，但是却以 减少变量的解析度为代价，说明初始因子建立的模型并不能解决问题，是不可取的，因 而仍需继续考虑旋转后的因子情况，在确保 全部的因子都可解析的条件下，进一步降维以找出影响北京市水资源短缺的主要风险因子。 第 9 页 /共 27 页 表 3转后因子载荷阵 变量 因子 1 因子 2 因子 3 因子 4 因子 5 旋转后因子载 荷阵计算编程看附录 尽管 的相关阵特征值的方差贡献率建议选 2 个或者 3 个初始因子进行建模，但是从旋转后的因子载荷阵可以看到，在第 4号和第 5号因子中仍然存在特别显著的特征值，因此并不能舍弃第 4号和第 5号因子，还需要 借助 因子载荷 阵 绝对值 0、 1 两极分化对比表 进行进一步的判断 。 定 主要风险 因子个数 子载荷阵绝对值 0、 1 两极分化频数对比表 借助 0因子载荷阵绝对值 0、 1 两极分化频数对比表 (如表 3表 3 因子载荷 阵 绝对值 0、 1 两极分化对比表 因子载荷区间 上 合计 频数 初始 5 4 5 1 65 80 旋转后 8 6 2 0 64 80 查找检验相关系数的临界值表， (r、 (r， 表 3子 第 10 页 /共 27 页 载荷阵绝对值 0、 1 两极分化对比表很明显的看到， 在初始因子载荷阵中，显著性大于 4个，而在旋转后的因子载荷阵中却有 16个；此外，在各段的显著水平的变量个数呈下降趋势。因此， 无论是从特别显著性还是从可解释变量来看，旋转后的因子载荷阵都优于初始因子载荷阵。 定 风险 因子个数 在表 3，显著性在 上的因子个数，旋转后的因子载荷阵比初始因子载荷阵多 3 个，在旋转后因子载荷阵的 5 个因子中，显著性在 上的可解释变量为 16 个，恰好是我们所选变量的个数，所以选择旋转后的因子作为我们的因子分析优化模型的 主要 因子 。 将 前 5 个旋转因子设为1f，2f，3f，4f，5f，在 变量正态分布下，取显著水平为 5%，由 7(子1f，2f，3f，4f，5f与变量显著相关；其它因子与变量没有显著相关 。因此我们选用旋转后的因子建立因子分析优化模型，并确定 风险 因子个数5m 。 险因子 命名 及其正向化 由 显著相关的临界值为 (子1f与7X、9X、11X、12X、13X、1510X、147X、11X、13X、15著性在 90%以上， 故称1f为环保效益因子，正向化1f取正号；因子2f与2X、3X、4X、6X、8与3X、6X、8著性在 90%以上 ，故称2f为自然需求因子， 正向化2f取正号。 因子3f与5与5著性在 90%以上 ， 故称3f为 供水 因子，正向化3f取正号。因子4f与16与16著性在 80%以上 ， 故称4f为万元 耗下降率因子，正向化4f取正号。 因子5f与1与10%以上 ， 故称5f气象因子，正向化5f取正号 。 第 11 页 /共 27 页 要风险 因子得分函数 考虑到 n p 时， 件过程命令此时计算因子的变量系数矩阵是错的，须用 块另外计算小样本解因子的系数矩阵，具体参考 文献 6。 通过对小样本解的计算，得到以下五个 主要风险 因子的得分函数参数估计表(如表 3 表 3转后的因子参数估计表 变量 因子 1 因子 2 因子 3 因子 4 因子 5 旋转后因子参数计算编程看附录 因此， 根据表 3到 五个 主要风险 因子 的 得分函数 (参数保留 3位有效数字 )分别是： 第 12 页 /共 27 页 合 风险 因子得分函数 在五个主要风险因子得分的基础上，鉴于 因子 间 互 不相关，综合起来可反映样品的因子累加综合状况 (不是反映多变量信息最大化时的样品值状况 )， 同时也可以得到综合风险评价的得分。因此， 以旋转因子信息贡献率 /i p为权数构造水资源短缺 综合 风险评价 因子函数 (参数估计如表 3 表 3合 风险评价 因子参数估计表 (参数 1 2 3 4 5 6 7 8 综合因子 数 9 10 11 12 13 14 15 16 综合因子 综合风险评价因子参数估计计算编程看附录 0 8 2 2 0 2 5 5 1 7 8 0 0 2 8 9 2 1 7 16/1 8 9 0 6 4 4 5 3 7 5 0 8 1 5 3 8 2 3 3 2 综 第 13 页 /共 27 页 旋转后 各 主要风险 因子 及 综合因子样品值排序 ， 见表 3 3转后因子得分、综合 风险 因子样品值及排序 年份 1f 排序 2f 排序 3f 排序 4f 排序 5f 排序 综2008 2009 2007 2006 2005 2004 2003 2002 2001 1 3 4 5 6 7 8 9 5 3 2 4 6 7 9 8 7 8 4 6 5 3 2 1 2 4 3 7 5 9 1 6 3 4 9 6 2 5 7 8 2 3 4 5 6 7 8 9 （旋转后各因子得分排序编程看附录 类分析图表 调用 件的聚类分析类 最短距离法 过程命令， 得到聚类过程表 3按 主要风险 因子 样品值 1f 、 2f 、 3f、 4f 、 5 表 3短距离法 聚类过程 （ T i e 8 2 . . . 2 . . . 2 . . . 4 . . 3 . . 7 . . 8 . . 9 . 旋转后各因子得分排序编程看附录 第 14 页 /共 27 页 0 20 30 40 50 60 70 类结果 通过 主要风险 因子 样品值 对 9 个样品进行聚类。 结合最短距 离法的聚类过程， 经分析 、验证 以及以往的经验 ， 我们一致认为 取阀值为 38 的分类最合适，因此 ， 可分为五类 如下 ： 第一类： 2001 第二类： 2002,2003 第三类： 2005 第四类： 2006,2007,2008,2009 第五类： 2004 类 结果综述 结合前 5个旋转后因子得分样品值的聚类分析结果，因子得分、综合因子得分样品值和排序，因 子得分、综合因子得分函数，原始数据，原始变量名称的意义，进行优势、劣势和影响因素等的综合评价， 我们得出结论如下： 根据 分类结果可知， 从 第一类样品 到第五类样品，它们的综合风险因子评价有递增的趋势 ，说明 水资源短缺 风险在逐渐降低 ， 也间接说明了人们已经逐渐意识到了水资源短缺的严重性，也在潜移默化的改变着自己的生活习惯，形成了自主节水的意识，与此同时，政府部门也抓紧了对水资源的调控 。 以第一类样品 2001年为例，其综合因子为 排名第 9，说明在这 9年里其水资源短缺 现象并不是非常 严重。1f排名第图 3类分析图 第 15 页 /共 27 页 2009 年地下水供水62%地表水供水13%再生水18%南水北调水7%2008年南水北调水2%再生水17%地表水供水16%地下水供水65%图片来源： 2009年北京市水资源公报 9，2f子 值为 排名第 8，3f子 值为 排名第 1，4f耗下降率 因子 值为 名第 9，5f名第 8。原因是1X、2X、3X、4X、5X、6X、7X、8X、9X、10X、11X、12X、13X、14X、15X、16、 3、 6、 7、 1、 6、 9、 6、8、 1、 9、 9、 9、 1、 9、 3，影响了 综合 因子的排名。 说明在 2001 年期间，人们的节水意识相当低，尽管水资源供应量 非常 充足 ，但是万元 耗下降率却从未降低 和水资源的污染加剧 ，以致经济发展 依 靠环境的 破坏为代价，从而导致了水资源短缺风险的存在。 我们 建议： 2001 年是 水资源短缺相对严重的年份 ，在供水和自然条件方面虽然方面有较大的优势，然而环保效益和 万元 有进一步的提升空间。只要继续保持自己原有的优势，并且致力于提高水资源利用效率，降低消耗， 充分发挥人工增雨的作用 ， 就能 降低水资源短缺风险。 以第 四 类样品 2006年至 2009年 中的 2008 年 为例，其综合因子为 名第 1，说明在这 9年里其水资源短缺 最严重 。1f名第 2，2f子 值为 排名第 1，3f子 排名第 9，4f子 值为 名第 8，5f名第 1。 主要 原因 是 在奥运年里，大量的人口涌入北京， 自然需求突然增大，大大超过了水资源供应的能力， 以致水资源 严重 供不应求 (从供水因子为负值 以及 2008 年、 2009 年供 、用 水量对比图都 可 以 看出 )。 我们 建议： 2008 年是奥运年，秉承绿色奥运的理念，在节水方面做得 相当出色 ， 在水资源利用和降雨量方面虽然有较大的优势，然而 供水 量 却严重不足，导致水资源短缺 最 严重 。因此，在 致力于 改善环境的同时，也需要 保持自己原有的优势， 同时 加大水利工程的建设 和建图 308、 09 年供水量对比 第 16 页 /共 27 页 2008年环境用水9%生活用水42%工业用水15%农业用水34%图片来源： 2009年北京市水资源公报 立水资源应急响应措施，以 提高 临时的 水资源供应， 这才 能 降低水资源短缺风险。 四、 预 测 和检验 次指数平滑法 应用此模型对 2010年北京市水资源的利用情况进行预测，由于 2010年的各个变量值 均 未知，故 需对1 16行估计。 考虑到一次指数平滑法既不需要存储的全部历史数据，也不需要存储一组数据，从而可以大大减少数据存储问题， 因此， 在此 采用一次指数平滑法对各变量的值进行估计。 一次指数平滑法是指以最后的一个第一次指数平滑 ，即根据前期的实测值和预测值，以加权因子为权数进行加权平均预测后期的值 。 它是用一个 比例系数来表示指数平滑值反映数据 序列的变化 ， 如果为了使指数平滑值敏感地反映最新观察值的变化，应取较大 值，如果所求指数平滑值是用来代表该时间序列的长期趋势值，则应取较小 值。 其基本的 计算 公式为 2, )1( 112 其中量的 预测值，量的 真实 ， 0,1 。 量值估计 考虑 预测不可避免具有误差， 采用一次指数平滑法时的一般 做法是取几个 求 得相应 的预测值 ， 根据求得的残差平法和 来判定取哪个 ，残差平方 和小的那个 即为所取值 ，相应求得的值即为预测值 。 基于一次指数平滑法，对 2010年各变量的值进行预测，分别 取 为 由残 差平方和最小原则可 作出如下 判定 ： 对变量1 12计 2009年环境用水10%生活用水41%工业用水15%农业用水34%图 308、 09 年 用 水量对比 第 17 页 /共 27 页 取 ， 而对13 16计 ，最终得到 估计 结果如表4相关程序见附录。由于篇幅关系，在此只给出 2010 年的预测值及前一年2009年 的真实值、 估计 值 。 ） 型 综合 风险 评价 预测 将各个变量 2010 年的估计值代入 的综合因子评价函数中，得到结果是 ，与 2009 年的相比， 2010 年的综合评价值有所下降， 说明政府 还在保持着奥运年的管理理念 ，但仍有很大的降低风险空间，同时也说明政府在奥运过后对北京市环境 (特别是水资源环境 )的治理有所放松，以致于 奥运年一年后在大