三峡水库的滞温评估

　　4.2三峡水库的滞温评估。

　　宜昌位于长江上游与中下游的过渡地区，中华鲟的产卵场也位于此。宜昌水文站检测断面处于三峡大坝-葛洲坝梯级水库的下游，本研究拟以宜昌江段为例，以宜昌江段水温过程为研究对象，发掘三峡水库滞温的时间效应，对三峡水库的滞温效应进行定量评估。

　　4.2.1 时间序列相似性匹配技术。

　　时间序列相似性问题是时间序列数据挖掘研究中比较基础性的问题[105].数据挖掘是指从大量数据中提取有效的、新颖的的、潜在有用的、最终可被理解的模式的非平凡过程[106].简单来说，数据挖掘就是从数据库中有目的性的提取有用的信息。数据挖掘是一个多学科领域，它融合了数据库技术、人工智能、机器学习、统计学、知识工程、信息检索等各种手段、方法，综合性和应用性都非常好。相似性问题的研究是进行相似性搜索、聚类[107]、异常点检测、模式发现[108]、分类分段[109]等重大技术的前提。相似性匹配技术应用非常广泛，如水文数据分析[110]、图像匹配[111]、信息检索[112]、过程监控[113]、故障排查[114]、地震预测[115]、脑电图分析[116]、生物信息比对[117]等。

　　生态水文学是一个数据密集型研究领域，长期的观测和时间积累了大量水文数据，这些数据中蕴含着自然界水文、生态系统、生物等长期的演变规律和人类活动影响的信息。利用人工智能和数据挖掘技术从这些大量数据中发现水文规律和生态学机理，是信息技术服务于其他研究领域的重要应用。目前水文时间序列数据挖掘技术正获得越来越多的关注和应用，但将数据挖掘用于生物与非生物环境关系的研究尚显不足，需更多的学者发现并关注。

　　时间序列相似性匹配从序列长度的角度可以分为两类[118]:全序列匹配和子序列匹配，也称基于整体的和基于局部的相似性匹配。全序列匹配特点是比较序列和被比较序列长度相等；子序列匹配则是被比较序列长度大于比较序列，目的是在长序列中查找包含与短序列相似的子序列。根据查询及结果的要求，相似性搜索又可以分为三种：（1）。 指定查询序列搜索：给定一个查询序列，在目标序列中查找所所有与之相似的序列或子序列；（2）。 匹配序列对搜索：在目标序列中发现所有彼此相似的序列对或子序列对；（3）。 最近邻搜索：给定一个查询序列，在目标序列中查找 K 个与之最相似的序列或子序列。

　　相似性匹配需要解决的两个主要问题是模式（特征）表示和相似性度量[119].

　　模式表示是对时间序列进行抽象和概括的特征表示方法。由于时间序列相似性匹配技术是随着计算机技术的发展而产生，大量的数据存储及空间的占用显然是不符合时代特性的。模式表示的的关键作用就是解决这个问题，其通过特诊抽取或转化，在起到降维的同时又能反映原始时间序列必要的信息。因此，基于特征的近似匹配实用性往往比精确匹配强，更受欢迎。

　　由于时间序列的相似性匹配是建立在时间序列相互间的距离上产生的技术手段，因此其核心问题是相似性度量标准。通常，时间序列模式表示方法都伴随着相应的时间序列相似性度量方法，用来度量时间序列间的相似性，且目前的主要度量都是基于时间序列特征的。时间序列的特征可以分为三种：形态特征，结构特征以及模型特征，也可以分为整体特征和局部特征[120].形态特征通常指时间序列的形状变化特征，包括整体特征和局部特征，如时间序列的趋势、转折等；结构特征是对时间序列全局构造或内在变化机制的描述，通常反映的是整体特征，如均值、方差等；模型特征是通过运用模型分析得到时间序列的一些参数特征，能够反映时间序列某些潜在的规律，包含整体特征和局部特征，如自相关性、混沌特性等。形态特征往往反映时间序列的特征比较直观，包含的信息较少；结构特征则能够一定程度上反映时间序列的内在结构，有利于进一步分析事物发展的规律；模型特征则能发现形态和结构特征无法表达的时间序列的一些特性，但这些特征往往缺乏事物发展的自然或科学含义。

　　常用的几种基于结构特征的度量方法有：欧式距离，动态时间弯曲（DWT）[121],Pearson 相关系数，Spearman 相关系数[122]等。下面给出这个度量方法的距离计算公式：（1）。 欧氏距离（4.1）（2）。 Pearson 相关系数距离（4.2）（3）。 DWT 距离（4.3）其中，类似。

　　（4）。 Spearman 相关系数距离（4.4）表示在时间序列 X 中的秩，类似。

　　不同的度量方法适用不同的类型，有不同的特性。例如，欧式距离对序列的相似性要求最严格，无法消除噪声对相似性的影响；Pearson 相关系数能够容许时间序列在尺度位移或拉伸状态下的相似，但不能适用强度位移或拉升的情形。宋辞等对上述常用的几种方法的相似性度量特性做了概括120.

　　随着时间序列相似性匹配技术的应用越来越广泛，许多新的改进方法正不断被发现，没有哪一种方法是最优的，使用者在选择相似性度量标准时应根据目标序列本身特性及数据分析所要达到的不同目的来选择合适的模式表示方法和度量标准，而且积累一定的先验知识往往比复杂的算法获得更好的匹配效果。

　　4.2.2 基于时间序列平移的分解模型。

　　周期性时间序列一般分为周期性成分和非周期性成分，非周期性成分通常又包含趋势成分和噪声成分。标准周期性时间序列（如图 4.1 中序列 a）不包含趋势成分，周期性成分可以表示为一个固定周期的不断重复（序列 b），噪声部分为白噪声序列（序列 c）。理想状态下自然界非生物信息都跟随太阳公转和地球自转形成的能量周期呈现这种发展模式，但自然条件下，由于受到生物非稳定性因素的干扰，自身的环境不再围绕外界能量一个因素运转，人们观测到的指标通常不再是这种平稳周期性时间序列，而是发展了的非标准周期性时间序列。

　　非标准周期性序列即一般周期性序列，包含周期成分、噪声成分和趋势成分三者。若将噪声成分限制为白噪声序列，即只代表随机波动成分，周期性成分和趋势成分则分别反映了结构变化和趋势变化。

　　从大的时间尺度上讲，一个周期性时间序列可以看成一个维度等于周期长度的多元时间序列，即每个周期序列单独作为多元时间序列在该时期的状态，多个按时间先后顺序的状态共同决定了周期性时间序列的结构和趋势。因此，研究周期性时间序列的周期间差异，是揭示时间序列发展的有效途径。

　　两个时间序列间的差异不仅是统计值（如均值）的差异，而且结构也可能存在差异。结构相同的序列间差异可以用其统计值的差值来衡量，即趋势变化；但结构变化的时间序列，由于变化形式众多，用差值甚至少数几个统计值很难快速、准确且全面反映序列的变化。

　　周期性时间序列平移分解模型的原理类似正弦曲线的平移理论，其适用于均值变化的趋势改变和相位变化的结构改变类型，两种变化都表现为在不考虑噪声影响情况下过程基本保持一致，仅在二维坐标系内发生整体平移。趋势变化表示周期内时间序列的每个值大小同步发生变化，即沿坐标系纵轴发生平移；相位变化表示周期过程沿坐标系横轴发生平移，如图 6.2-2 中所示，d 序列为对照序列，序列 e 为趋势变化情形，序列 f 为相位变化情形，周期性时间序列平移分解模型适用于每个周期内只存在这两种情形以及这两种变化同时存在的情形。

　　从整个周期性时间序列来说，周期性表现为时间序列在每个周期内的过程都保持固定的或者相似的模式，每个周期都是在遵从这个固定模式的基础上在对应时间位置，即相隔一个或若干个周期长度，发生一定的形变和波动，从而形成趋势甚至结构变化。当这些变化可以用上述的趋势变化和相位变化来描述时，那么整个时间序列的所有的周期变化都可以用两个参数来衡量：趋势位移和相位移，这个两个参数构成的序列则可以简单、有效地反映该周期性时间序列的趋势和结构变化过程。

　　因此，周期性时间序列平移分解模型主要分为两个步骤，一是固有周期模式的确定；二是计算趋势位移和相位移。下面针对时间序列周期过程的变化，分三种类型介绍：

　　（1）。 趋势位移；发生趋势位移的周期性时间序列，固有周期模式只在周期内同步增减，因此取多个周期的平均过程作为固有周期模式，便能反映该周期性时间序列的发展模式。发生趋势位移的周期过程对应的首末端点位置并未发生移动，因此位移距离可以为任意对应点的差值，为消除噪声影响，计算时采用周期过程和国有周期模式的均值之差。

　　（2）。 相位移：

　　发生相位移的周期性时间序列，平移后的周期过程出现错位，由于噪声影响，平移后序列结构缺乏明显标志，周期的起止时间点也变得模糊，使得平移后的位置难以捕捉，也难以确定固有周期模式。为综合考虑整个时间序列的周期过程，考虑选取研究期内与所有周期最相似的周期过程作为国有周期过程。

　　相似性度量采用欧式距离平方，如公式 9,固有周期过程计算两两周期过程间相似度，最后选定所有周期中与其他周期过程相似性度最高的作为固有周期过程。

　　再利用周期过程与固有周期模式的相似性，通过最大相似度捕捉到平移后的固有周期序列，再计算每个周期对应的平移距离。此时计算的相似度是序列中所有的长度等于一个周期的子序列与固有周期模式的相似度，每个周期都对应一个最大相似度，其对应的子序列即作为平移后的固有周期序列，与实际周期过程的位移距离即为相位移。

　　（3）。 上述两种位移同时存在：

　　当同时存在趋势位移和相位移时，分解方法同上述第二种情形，但在计算相似度前，应先对目标序列进行标准化，标准化方法如下：

　　（4.5）式中为目标序列的均值。标准化可以首先忽略趋势位移，在确定了相位移距离及平移后的固有周期位置以后，再利用两序列均值差确定趋势位移距离。

　　序列长度较长时，应先检测周期时间序列的结构突变点，不同的周期结构发生的结构变化类型可能不同，发生上述两种情势的位移变化时才能采用该分解模型，并以最初阶段的固有周期模式作为整个序列的固有周期模式再进行下一步分析。

　　4.2.3 葛洲坝下周期水温过程的分解。

　　4.2.3.1 模型适用性分析影响葛洲坝下游宜昌江段水温过程变化的主要因素为气候变化及长江上游梯级水库的运行。1981 年-2012 年水温过程不仅发生了趋势变化，周期也发生了变化。

　　通过突变分析，确定水温过程突变发生在 2005 年至 2006 年间，且主要突变类型为升温过程的推迟及降温过程的滞后，初步推测造成这种突变的主要原因是三峡水库的蓄水。为定量揭示气候变化和长江上游梯级水库的运行分别对长江中游水稳过程的影响，引入时间序列平移分解模型。

　　不妨猜想：

　　a. 葛洲坝下游宜昌江段水温过程的固有周期模式不发生结构变形。

　　b. 气候变化仅影响固有周期模式的纵向平移。

　　c. 水利工程建设仅影响固有周期模式的横向平移。

　　现给出证明：

　　（1）。 猜想 a:

　　不难发现，由于水温过程在去除气候变化和水利工程建设两者主要影响因素外，也就是说，猜想 2 和 3 同时成立的条件下，猜想 1 就自然而然的成立。

　　（2）。 猜想 b:

　　气温与水温关系早有研究证明，未受人类干扰的自然河流，水温与同气温同步变化，仅滞后 2h 相位。因此，去除水利工程建设的影响，气候变化对水温过程的固有周期的影响基本为线性相应关系，若气候过程为平稳周期性时间序列，则猜想成立。

　　分析选用宜昌站 1981-2012 年 32a 逐月平均气温过程进行周期平稳性验证，方法如下：

　　首先进行去周期性操作。计算公式为：

　　（4.6）对残差序列进行 mann-kendall 趋势检验，显著性水平选 5%.结果表明残差序列存在上升趋势，说明该研究期内气温过程发生了趋势变化。

　　再对残差序列进行一阶差分，消除趋势性成分，继续对残差序列进行趋势检验。

　　结果表明残差序列为平稳时间序列，不存在趋势性，说明该时期内气温过程只存在趋势变化，不存在结构变化。因此，猜想 b 成立。需要注意的是：尽管纵向位移全部由气候因素影响产生，但不等同于气候变化效应。应为气温对水温的影响本身并没有滞后。

　　（3）。 猜想 c在第二章节中，我们已对葛洲坝下游水温过程进行结构变异诊断，发现 2005年至 2006 年间出现了突变，突变类型变现为升温过程的推迟及降温过程的滞后。

　　若强行以周期模式发生改变来描述该水温过程的变化，确有不妥，但仔细观察发现 2006 年以后，这种突变导致了一年 12 个月中 8-9 个月水温过程滞后，最低温月份推迟到了 3 月，只有高温月份未显著滞后。且研究期内水温过程的结构突变主要由水利工程建设这一因素引起，这使得用固有水温周期过程的横向平移来评估滞温效应的可行性大大提高。该分解过程不能完全覆盖水温所发生的结构变化，但能够为评估水利工程建设对水温过程的影响提供较直观的评价。因此，猜想 c在此种目的意义下也可以成立。

　　猜想的全部成立，充分保证了利用时间序列平移分解模型进行两种因素对水温过程的综合影响结果进行分解的前提条件，即时间序列的固有周期模式不发生结构变化。

　　4.2.3.2 分解过程。

　　计算工具为 Matlab2014b,计算过程如下：

　　1. 计算 1981 年~2005 年多年平均逐日水温，得到长度为 365 的对照序列 T0,如图 4.4;根据上一节对 32 年水温过程的结构诊断，突变发生在 2005 年后，因此初步认定 2005 年前，水温过程基本由气候变化影响，因此固有周期过程必须从 1981 年~2005 年内提取。

　　由于气候因素只产生趋势变化，因此认定，1981 年-2005 年内水温过程也只产生了趋势变化，未发生结构变化，可根据第一种情形选定固有周期模式。本研究采用 1981 年~2005 年多年平均逐日水温过程作为固有周期过程。由于这 25 年内，气温主要成单调升高趋势，因此在距离度量中，不会损失气候变化的影响结果，而且不会影响横向平移的分解。

　　2. 以欧式距离，分别计算 1982 年~2011 年水温过程与对照序列的最大相似度。

　　以 1982 年为例，以 1 月 1 日为参照起点，运用公式 4.1,依次计算水温序列：

　　T（1）：1982/1/1~1982/12/31,T（2）：1982/1/2~1983/1/1,T（3）：1982/1/3~1983/1/2,……与对照序列 T0的相似度 S（1），S（2），S（3）……；相似度计算前对两条序列分别作标准化处理，方法见公式 4.5.

　　同时计算序列：

　　T（-1）：1981/12/31~1982/12/30,T（-2）：1981/12/30~1982/12/29,T（-3）：1981/12/29~1982/12/28,……与对照序列 T0的相似度 S（-1），S（-2），S（-3）……，标准化过程同上。

　　一般可初步推至 T（40）和 T（-40），得到长度为 80 的相似度序列 S,求相似度序列 S 的极小值 min（S）。

　　3. 分别求 1982 年~2011 年水温过程的滞温时间。

　　求得相似度序列以后，相似度极小值 min（S）对应的序列的实际日期与 1 月 1 日的差值即为滞温时间。例如，1982 年的相似度最小值如果对应的日期为1982/1/3~1983/1/2,则其滞温时间为 2 天，代表其周期水温过程由 1 月 1 日推迟至1 月 3 日；若相似度最小值对应的日期为 1981/12/29~1982/12/28,则其滞温时间为-2 天，代表其周期水温过程由 1 月 1 日提前到了 1981 年的 12 月 29 日。如果极小值出现在 T（40）或 T（-40），可将计算范围拓展至 T（50）和 T（-50）或者更长。

　　4.2.3.3 分解结果。

　　2003 年以前，滞温时间相对稳定且时间不长，基本不超过 10 天，可以认为这一段时期宜昌水温过程没有发生滞后现象，2003 年以后，宜昌江段水温滞后时间有了明显增长，且随着年份增加，滞温时间也随之增加。为提高模型准确性，将固有周期水温过程调整为 1981 年-2003 年 23a 平均逐日水温过程，模拟结果与 25a 结果基本一致，原因是 2003 年-2005 年滞温效应并不显著，这也是水温过程突变发生 2005 年的主要原因。

　　模型结果与模型选择的对照期非常吻合，表明选定 1982 年~2003 年作为对照期是合理且正确的，模拟结果最终说明：长江上游梯级水库运行对长江中游地区的水温过程影响发生在 2003 左右，滞温效果非常显著，并且滞后时间持续增加，到 2011 年滞后期已接近 30 天。

　　图 4.6 是模拟的 2011 年的水温滞后效果图，可以看到，2011 年葛洲坝宜昌江段水温相对于三峡大坝蓄水前有了较大幅度的滞后，通过基于时间序列平移的分解模型，对水温过程的升温、降温以及最高温、最低温都有较好的模拟。

　　4.2.4 分析与讨论。

　　三峡大坝位于葛洲坝上游 38 公里，2003 年启动下闸蓄水，2006 年水位提升至 156m,2010 年三峡水库试验性蓄水首次达到 175m 最终水位，总库容量达到 393亿立方米。三峡水库的水位变化，将显著影响水库的水体层化和水温状态，同时改变水气界面的热交换，这些都会影响三峡大坝的下泄水温。三峡水库属于弱分层型水库，但水库的结构及运行状态对水库水温起关键调控作用，宜昌江段中华鲟产卵场水温也主要受其影响。从 2003 年开始，三峡水库的蓄水就对长江中下游水温过程产生了影响，但水位还处于较低水平，水库的表底水温差异较小，水温过程的滞后也相对较短。从 2003 年到 2011 年，三峡水利枢纽还处于试运行期，水位在逐步提高，滞温效应也逐渐增长，水温过程滞后期的延长是三峡大坝建成蓄水后的累积效益。但就目前来说，该趋势将是短期的，随着三峡水利枢纽逐渐从试运行过渡期转向稳定运行期，三峡水库的温度模式也将逐渐稳定，下游水温的滞后也将随之稳定。

　　三峡水库的滞温效应将对长江流域中下游的所有水生生物的栖息和繁衍造成巨大的影响。生活于长江流域的大多数鱼类对水温有原本的适应性，生殖活动也一般发生在固定的季节。鱼类原产卵场的水温过程不断推迟，本来鱼类的繁殖随之活动向后推迟也非常简单，但同时鱼类繁殖所需的其他水文条件并未推迟，这就造成了鱼类繁殖行为将重新对这种不同的水环境条件进行适应，这种节律性的改变是对生物自身调节能力的最大考验，适应能力弱的生物无法及时调整应对，种群数量必然会在这场变迁中很快衰减，进而流域的生物群落结构也将发生改变，生物的多样性及生物量短期内下降是必然的，而一些关键物种能否保存下来往往是焦点问题。

　　4.3小结。

　　本章介绍了大型水库的水温层化理论、水库对下游水温的影响机理以及三峡水库滞温效应的研究现状，指出三峡水库的滞温是影响长江流域水生态和水生生物的重要环境因素，目前正受到越来越多关注，但定量研究相对欠缺。

　　结合时间序列相似性匹配技术，本章建立了基于时间序列平移的分解模型，能够适用于发生模式平移的结构变化的周期性时间序列，识别周期性时间序列中模式偏移现象，并准确量化模式的纵向和横向平移距离。利用该模型，对三峡水库的滞温效应（仅宜昌江段）做了定量模拟，模拟效果与模型假设即数据水文背景向吻合，不仅找出了滞温产生的时间起点，而且清晰显示了建坝以来三峡水库产生的年水温过程滞后时间的变化趋势。