CN112035553B - 一种基于空间统计理论的水文设计值计算方法 - Google Patents

一种基于空间统计理论的水文设计值计算方法 Download PDF

Info

Publication number
CN112035553B
CN112035553B CN202010911383.3A CN202010911383A CN112035553B CN 112035553 B CN112035553 B CN 112035553B CN 202010911383 A CN202010911383 A CN 202010911383A CN 112035553 B CN112035553 B CN 112035553B
Authority
CN
China
Prior art keywords
tested
station
design
hydrological
calculating
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010911383.3A
Other languages
English (en)
Other versions
CN112035553A (zh
Inventor
余挺
夏传清
康有
马顺刚
张军良
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
PowerChina Chengdu Engineering Co Ltd
Original Assignee
PowerChina Chengdu Engineering Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by PowerChina Chengdu Engineering Co Ltd filed Critical PowerChina Chengdu Engineering Co Ltd
Priority to CN202010911383.3A priority Critical patent/CN112035553B/zh
Publication of CN112035553A publication Critical patent/CN112035553A/zh
Application granted granted Critical
Publication of CN112035553B publication Critical patent/CN112035553B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F16/00—Information retrieval; Database structures therefor; File system structures therefor
    • G06F16/20—Information retrieval; Database structures therefor; File system structures therefor of structured data, e.g. relational data
    • G06F16/24—Querying
    • G06F16/245—Query processing
    • G06F16/2458—Special types of queries, e.g. statistical queries, fuzzy queries or distributed queries
    • G06F16/2462—Approximate or statistical queries
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F16/00—Information retrieval; Database structures therefor; File system structures therefor
    • G06F16/20—Information retrieval; Database structures therefor; File system structures therefor of structured data, e.g. relational data
    • G06F16/29—Geographical information databases

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Databases & Information Systems (AREA)
  • Probability & Statistics with Applications (AREA)
  • General Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Fuzzy Systems (AREA)
  • Computational Linguistics (AREA)
  • Software Systems (AREA)
  • Mathematical Physics (AREA)
  • Remote Sensing (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明涉及水利水电工程中的水文分析技术领域，本发明旨在解决现有的缺资料地区水文设计值计算较为困难的问题，提出一种基于空间统计理论的水文设计值计算方法，技术方案概括为：选取关键属性因子集；计算已测站点指定设计频率下的第一水文设计值；选取回归法因子，建立区域回归方程，计算已测站点的第一地区综合系数以及回归方程参数；选取普通克里金法因子，建立基于属性因子的普通克里金法，计算待测站点的第二地区综合系数；计算待测站点指定设计频率下的第二水文设计值。本发明实现了缺资料地区的水文分析计算，并且其准确性和可靠性较高。

Description

一种基于空间统计理论的水文设计值计算方法
技术领域
本发明涉及水利水电工程中的水文分析技术领域，具体涉及一种基于空间统计理论的水文设计值计算方法。
背景技术
水文设计值准确可靠计算对于水利水电工程规划设计、施工建设和运行调度等至关重要，是水利水电工程设计所要面临和解决的首要任务。近年来，水利水电工程规划建设整体呈现有资料地区转向缺资料地区的发展态势。而当地水文基本资料往往极为匮乏，已成为限制当地水资源开发利用的技术瓶颈之一。资料匮乏地区工程水文设计所面临的困难越来越大，而水文设计值对合理确定工程规模及确保工程安全可靠具有极其重要的意义。传统工程水文分析方法在资料匮乏地区难以准确可靠计算水文设计值，其成果具有较大的不确定性，一直是长期困扰国际水文学界的难题之一。缺资料地区水文预测研究已成为国际水文及水资源领域研究的前沿和热点问题之一，具有重要的理论意义和实践价值。2003年，国际水文科学协会(International Association of Hydrological Sciences，IAHS)启动了第二个国际水文计划，确定了未来十年重大研究计划——“无资料地区水文学研究”(Prediction in Ungauged Basins，PUB)，为缺资料地区工程水文分析计算开拓了新的机遇及发展方向。
目前，国内外在缺资料地区水利水电工程水文设计值计算中常采用水文比拟法，首先选择与研究流域属性相似且具有长系列水文资料的邻近流域作为参证流域，然后将参证流域的水文信息(如方法参数值、水文设计值、水文特征值等)直接或经过修正后移植至研究流域。由此可见，水文比拟法一般侧重于寻找邻近相似流域，以流域间的相似性为基础，实现相似流域间水文信息的传递。因此该方法的前提在于研究流域的邻近地区存在属性相似且具有长系列水文资料的流域。但是，在实际中缺资料地区的邻近流域可能属性不相似，其水文过程差异较大，造成“相似的无资料、有资料不相似”的困境。因此，寻找具有长系列水文资料的邻近相似流域在实际中未必是一种可行的途径。如何基于高质量的多源海量遥感大数据，充分利用属性差异较大的邻近流域水文设计成果，是缺资料地区水文设计值计算亟待解决的关键科学问题，在水利水电工程规划设计中具有重要的应用价值。
发明内容
本发明旨在解决现有的缺资料地区的水文预测较为困难的问题，提出一种基于空间统计理论的水文设计值计算方法。
本发明解决上述技术问题所采用的技术方案是：一种基于空间统计理论的水文设计值计算方法，包括以下步骤：
步骤1、基于可获取的多源海量遥感大数据，提取研究流域及其邻近流域的多种属性栅格数据；建立基于地理探测器的属性因子优选方法，从所述属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集；
步骤2、收集并整理多个已测站点的水文资料，根据所述水文资料并建立基于SCE-UA算法的水文频率分析方法，计算各已测站点指定设计频率下的第一水文设计值；
步骤3、对于各已测站点，采用逐步回归法，从所述关键属性因子集中选取对第一水文设计值影响显著的属性因子作为回归法因子；根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，计算各已测站点的第一地区综合系数以及回归方程参数；
步骤4、对于各已测站点，采用逐步普通克里金法，从所述关键属性因子集中选取对第一地区综合系数影响显著的属性因子作为普通克里金法因子；建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数；
步骤5、根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，并基于所述区域回归方程计算待测站点指定设计频率下的第二水文设计值。
进一步的，还包括：
步骤6、建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数，分别计算各已测站点的第一综合系数的预测值；
步骤7、根据各已测站点的第一地区综合系数的预测值、对应的回归法因子以及回归方程参数，分别计算各已测站点的第一水文设计值的预测值；
步骤8、根据各已测站点的第一水文设计值及其预测值计算第一水文设计值的误差值，建立基于普通克里金法因子的反属性插值方法，根据所有已测站点的第一水文设计值的误差值计算待测站点的第二水文设计值的误差值；
步骤9、根据所述第二水文设计值的误差值计算待测站点的第二水文设计值的修正值。
进一步的，还包括：
步骤10、重复步骤1至9，计算待测站点所有设计频率下的第二水文设计值的修正值；
步骤11、建立基于SCE-UA算法的水文频率分析方法，计算待测站点的水文频率分布参数；
步骤12、由水文频率分布线型及其参数值，计算待测站点所有设计频率下的第二水文设计值的最终值。
进一步的，还包括：
步骤13、针对指定设计频率下，从所有已测站点中抽取一站作为验证站，剩余已测站点组合作为率定站；
步骤14、假设率定站和验证站分别为已测站点和待测站点，根据率定站指定设计频率下的第一水文设计值，计算验证站指定设计频率下的第二水文设计值的预测值；
步骤15、重复步骤13至14，将所有已测站点逐一抽取作为验证站，并将验证站假设为待测站点；依次取不同的设计频率，计算所有验证站所有设计频率下的第二水文设计值的预测值，并将其与对应设计频率下的第二水文设计值进行比较，进而评估该方法的可靠性。
进一步的，步骤8中，所述第二水文设计值的误差值
Figure BDA0002663415900000031
的计算公式为：
Figure BDA0002663415900000032
其中，
Figure BDA0002663415900000033
式中，
Figure BDA0002663415900000034
为第i个已测站点指定设计频率下的第一水文设计值的误差值，di为待测站点与第i个已测站点之间的属性距离，
Figure BDA0002663415900000035
为待测站点的属性坐标系，
Figure BDA0002663415900000036
Figure BDA0002663415900000037
为第i个已测站点的属性坐标系，B1,B2,…,Bk为普通克里金法因子，n为已测站点的个数，k为普通克里金法因子个数。
进一步的，步骤1中，所述从属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集具体包括：
采用地理探测器中的因子探测法，将水文变量栅格数据作为因变量，其余属性栅格数据作为自变量，探测因变量的空间分异性以及各自变量多大程度上解释了因变量的空间分异，选取对水文变量影响显著的属性组成关键属性因子集。
进一步的，步骤2中，所述计算各已测站点指定设计频率下的第一水文设计值具体包括：
假定水文变量X服从P-III型分布，记作X～Γ(x；a,α,β)，其概率密度函数f(x)和分布函数F(x)计算公式为：
Figure BDA0002663415900000038
式中，a为位置参数，α为形状参数，β为尺度参数；
根据指定的水文频率分布线型即P-III型分布，建立基于SCE-UA算法的水文频率分析方法，计算水文频率分布参数；
根据水文频率分布线型及其参数值，计算各已测站点指定设计频率下的第一水文设计值QP，其计算公式为：
QP＝Ex×(1+ΦP×Cv)；
其中，
Figure BDA0002663415900000041
式中，P为设计频率，Ex为期望值，Cv为变差系数，Cs为偏态系数，ΦP为离均系数，Γ为Gamma函数。
进一步的，步骤3中，所述区域回归方程采用多元幂函数方程，其计算公式为：
Figure BDA0002663415900000042
式中，QP为指定设计频率下的水文设计值，CP为地区综合系数，A1,A2,…,Am为回归法因子，α1,α2,…,αm为回归方程参数，m为回归法因子个数。
进一步的，步骤4中，所述建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数具体包括：
假设各已测站点的第一地区综合系数分别为CP(xi)(i＝1,2,…,n)，则待测站点的第二地区综合系数CP(x0)是所有已测站点的第一地区综合系数CP(xi)的线性加权和，即：
Figure BDA0002663415900000043
式中，CP(x0)为待测站点的第二地区综合系数，CP(xi)为第i个已测站点的第一地区综合系数，
Figure BDA0002663415900000044
为待测站点的属性坐标系，
Figure BDA0002663415900000045
为第i个已测站点的属性坐标系，B1,B2,…,Bk为普通克里金法因子，λi为第i个已测站点的克里金权重，n为已测站点的个数，k为普通克里金法因子个数。
进一步的，所述第i个已测站点的克里金权重λi的计算公式为：
Figure BDA0002663415900000046
式中，
Figure BDA0002663415900000047
为第i个和第j个已测站点xi与xj之间的半方差，
Figure BDA0002663415900000048
为第j个已测站点xi与待测站点x0之间的半方差，μ为拉格朗日乘子。
本发明的有益效果是：本发明提供的一种基于空间统计理论的水文设计值计算方法，充分利用多源海量遥感大数据，建立水文变量与流域属性因子之间空间定量关系，将有资料地区的水文设计值空间移植至无测站流域，实现缺资料地区的水文分析计算，有效提高了缺资料地区水文设计值计算结果的准确性和可靠性，具有可靠性高、普适性强、计算效率高的优势，为资料匮乏地区水利水电工程水文分析计算提供了一种更加稳健可靠的计算方法，为更好地进行缺资料地区水利水电工程规划设计提供了科学依据，具有较好的应用前景。
具体实施方式
本发明旨在解决资料匮乏地区水文设计值计算的难题，提供一种基于空间统计理论的水文设计值计算方法，包括以下步骤：步骤1、基于可获取的多源海量遥感大数据，提取研究流域及其邻近流域的多种属性栅格数据；建立基于地理探测器的属性因子优选方法，从所述属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集；步骤2、收集并整理多个已测站点的水文资料，根据所述水文资料并建立基于SCE-UA算法的水文频率分析方法，计算各已测站点指定设计频率下的第一水文设计值；步骤3、对于各已测站点，采用逐步回归法，从所述关键属性因子集中选取对第一水文设计值影响显著的属性因子作为回归法因子；根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，计算各已测站点的第一地区综合系数以及回归方程参数；步骤4、对于各已测站点，采用逐步普通克里金法，从所述关键属性因子集中选取对第一地区综合系数影响显著的属性因子作为普通克里金法因子；建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数；步骤5、根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，并基于所述区域回归方程计算待测站点指定设计频率下的第二水文设计值。
首先，提取研究流域及其邻近流域的多种属性栅格数据，从多种属性栅格数据中选取关键属性因子集；对多个已测站点的水文资料进行水文频率分析得到第一水文设计值，从关键属性因子集中选取对第一水文设计值影响显著的属性因子作为回归法因子；然后，根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，并采用最小二乘法计算回归方程参数；将各已测站点的第一水文设计值和对应的回归法因子代入区域回归方程中，计算得到各已测站点的第一地区综合系数；再然后，从关键属性因子集中选取对第一地区综合系数影响显著的属性因子作为普通克里金法因子，建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算用于对待测站点的第二水文设计值计算的区域回归方程中的第二地区综合系数；对于已测站点和待测站点来说，所述区域回归方程的回归方程形式和回归方程参数均相等；最后，根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，代入到所述区域回归方程中，即可计算得到待测站点指定设计频率下的第二水文设计值。
实施例
本实施例旨在解决资料匮乏地区水文设计值计算的难题，提供一种基于空间统计理论的水文设计值计算方法，其计算步骤包括：属性因子优选、水文频率分析、地区综合分析、空间插值分析、水文设计值预测、水文设计值修正、水文设计值优化、方法可靠性评估；其核心思想是：充分利用多源海量遥感大数据，建立水文变量与流域属性因子之间空间定量关系，将有资料地区的水文设计值空间移植至无测站流域；其实质是利用区域内多个已测站点的水文信息对待测站点进行空间估计；下面对实施例的各个步骤进行详细说明。
(一)属性因子优选
步骤1、基于可获取的多源海量遥感大数据，提取研究流域及其邻近流域的多种属性栅格数据；建立基于地理探测器的属性因子优选方法，从所述属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集；
所述多源海量遥感大数据是指国内外可获取的覆盖研究流域及其邻近流域边界范围的多种卫星遥感属性栅格数据，包括：径流、洪水、降水、蒸发、地形、土壤、植被等；
比较有代表性的卫星遥感属性栅格数据为美国国家航空航天局(NASA)发射的遥感卫星Terra和Aqua的中分辨率成像光谱仪(MODIS)生成的各种陆地标准产品，例如归一化植被指数(Normalized Difference Vegetation Index，NDVI)、白天地表温度(LandSurface Temperature Day，LSTD)、夜晚地表温度(Land Surface Temperature Night，LSTN)、叶面积指数(Leaf Area Index，LAI)、地表反射率(Surface Reflectance，SR)等，其最大空间分辨率可达250m。
所述提取研究流域及其邻近流域的多种属性栅格数据具体包括：
1)利用ArcGIS软件对卫星遥感属性栅格数据进行拼接、裁剪、投影和提取等操作，提取研究流域及其邻近流域的径流、洪水、降水、蒸发、地形、土壤、植被等属性栅格数据；
2)基于数字高程模型(Digital Elevation Model，DEM)栅格数据，利用ArcGIS软件分析提取研究流域及其邻近流域的经度、纬度、坡度、坡向、地形指数、地形开阔度、地表粗糙度、地形起伏度等属性栅格数据。
所述建立基于地理探测器的属性因子优选方法，从所述属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集具体包括：基于水文变量(指径流或洪水)栅格数据和降水、蒸发、地形、土壤、植被等其余属性栅格数据，采用地理探测器中的因子探测法，将水文变量栅格数据作为因变量，其余属性栅格数据作为自变量，探测因变量的空间分异性以及各自变量多大程度上解释了因变量的空间分异，选取对水文变量影响显著的属性组成关键属性因子集。
其中，地理探测器中因子探测法的基本原理是当被解释变量和影响因子的空间分异格局趋同，说明二者具备统计关联性，探测因变量的空间分异性以及各自变量多大程度上解释了因变量的空间分异；用q值度量，其值越大说明因变量的空间分异性越明显，其计算公式为：
Figure BDA0002663415900000071
其中，
Figure BDA0002663415900000072
SST＝Nσ2；
式中，h＝1,2,…,L为因变量Y或因子X的分层，即分类或分区，Nh和N分别为层h和全区的单元数，
Figure BDA0002663415900000073
和σ2分别是层h和全区的Y值的方差，SSW和SST分别为层内方差之和和全区总方差。
(二)水文频率分析
步骤2、收集并整理多个已测站点的水文资料，根据所述水文资料并建立基于SCE-UA算法的水文频率分析方法，计算各已测站点指定设计频率下的第一水文设计值；
其中，所述多个已测站点指位于研究流域及其邻近流域内至少收集到10个具有长系列水文资料的水文站；所述长系列水文资料指水文站的年平均流量、年最大洪峰流量、年最大时段洪量等资料，且系列长度应在30年以上。
所述建立基于SCE-UA算法的水文频率分析方法，计算各已测站点指定设计频率下的第一水文设计值，具体包括：
(1)结合我国水利水电工程水文分析计算实际情况，假定水文变量X服从P-III型分布，记作X～Γ(x；a,α,β)，其概率密度函数f(x)和分布函数F(x)计算公式为：
Figure BDA0002663415900000074
式中，a为位置参数，α为形状参数，β为尺度参数。
(2)根据上述指定的水文频率分布线型即P-III型分布，建立基于SCE-UA算法的水文频率分析方法，计算水文频率分布参数；其中，优化变量为Ex、Cv与Cs，目标函数为离差绝对值和，约束条件为参数取值范围，绘点位置采用现行《水利水电工程设计洪水计算规范》(SL 44-2006)中数学期望公式进行计算；
(3)根据水文频率分布线型及其参数值，计算指定设计频率的水文设计值，如指定设计频率的设计径流、设计洪峰或设计洪量等，其计算公式为：
QP＝Ex×(1+ΦP×Cv)；
其中，
Figure BDA0002663415900000081
式中，P为指定设计频率，Ex为期望值，Cv为变差系数，Cs为偏态系数，ΦP为离均系数，Γ为Gamma函数。
所述目标函数为离差绝对值和，即参数估计值所对应的理论频率曲线和绘点位置的纵坐标之差的绝对值和，其计算表达式为：
Figure BDA0002663415900000082
式中，S(θ)为目标函数值，θ为优化变量(包括Ex、Cv与Cs)，Qi为从大到小排列后第i个水文系列实测值，F-1(Pi；θ))为给定设计频率Pi所对应的水文设计值，即理论频率曲线的纵坐标值，L为水文系列样本容量。
所述约束条件为参数取值范围指采用线性矩法初估三个统计参数Ex、Cv与Cs的初始值Ex0、Cv0与Cs0，确定参数取值范围应该满足如下关系：
g1：Ex∈[Ex0×0.3,Ex0×1.7]；
g2：Cv∈[Cv0×0.4,Cv0×1.6]；
g3：Cs∈[Cs0×0.5,Cs0×1.5]；
g4：Cs∈[0,2]；
g5：
Figure BDA0002663415900000083
式中，Ex0、Cv0、Cs0为采用线性矩法初估的期望值、变差系数、偏态系数，xmin为实测水文系列中的最小值。
所述SCE-UA(Shuffled Complex Evolution University of Arizona)算法结合了随机搜索、传统复合型法和生物竞争优胜劣汰等方法的优点，是一种可以有效解决复杂的高维参数、多极值等非线性全局优化问题的进化算法，具有计算效率高、收敛速度快和求解稳定性好等优点，广泛应用于非线性复杂的水文预测方法参数率定中。因此，本发明中优选采用SCE-UA算法对三个优化变量(包括Ex、Cv与Cs)进行参数求解。
进一步地，在步骤(2)中，参数优化求解过程具体包括：首先将实测水文系列按照从大到小的顺序排列成x(1)≥x(2)≥…≥x(L)；其次根据现行《水利水电工程设计洪水计算规范》(SL 44-2006)中数学期望公式，计算各个样本x(i)对应的经验频率P(i)；然后在几率格纸上描绘出各个样本点据(P(i),x(i))，即为绘点位置；最后根据SCE-UA算法和目标函数优选出一条拟合效果最佳的曲线，其对应的参数值即为理论频率分布的参数最优值。
(三)地区综合分析
步骤3、对于各已测站点，采用逐步回归法，从所述关键属性因子集中选取对第一水文设计值影响显著的属性因子作为回归法因子；根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，计算各已测站点的第一地区综合系数以及回归方程参数；
可以理解，对于各已测站点，采用逐步回归法，按照流域属性因子对已测站点指定设计频率下的第一水文设计值的影响显著程度，从大到小依次逐个引入回归方程中，当先选入的属性因子由于后面的属性因子引入而失去重要性，就将它剔除；并且每引入或剔除一个属性因子都要作相应的F检验，直至方程既不能引入也不能剔除，即F检验通不过为止，则最终的属性因子作为回归法因子；
具体而言，区域回归方程采用多元幂函数方程，其计算公式为：
Figure BDA0002663415900000091
式中，QP为指定设计频率下的水文设计值，CP为地区综合系数，A1,A2,…,Am为回归法因子，α1,α2,…,αm为回归方程参数，m为回归法因子个数。
根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，并采用最小二乘法计算回归方程参数；将各已测站点的第一水文设计值和对应的回归法因子代入区域回归方程中，计算得到各已测站点的第一地区综合系数，其计算公式为：
Figure BDA0002663415900000092
(四)空间插值分析
步骤4、对于各已测站点，采用逐步普通克里金法，从所述关键属性因子集中选取对第一地区综合系数影响显著的属性因子作为普通克里金法因子；建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数；
具体而言，对于各已测站点，采用普通克里金法，按照流域属性因子对已测站点的第一地区综合系数的影响显著程度，从大到小依次逐个引入克里金方程中，当先选入的属性因子由于后面的属性因子引入而失去重要性，就将它剔除；并且每引入或剔除一个属性因子都要作相应的精度评价，直至方程既不能引入也不能剔除，即方法精度最高，则最终的属性因子作为普通克里金法因子。
普通克里金法是以变异函数理论和结构分析为基础，在有限区域内对待测站点的区域化变量进行线性无偏最优估计的一种优良方法，其估计方差最小。本发明建立基于属性因子的普通克里金法，将普通克里金法中坐标系采用优选出的普通克里金法因子，实现了将流域属性因子纳入克里金权重估计中，赋予流域属性相近的已测站点更高的权重，以估计待测站点的地区综合系数，进而可明显提高水文设计值计算结果的可靠性。
第二地区综合系数的计算方法如下：
假设各已测站点的第一地区综合系数分别为CP(xi)(i＝1,2,…,n)，则待测站点的第二地区综合系数CP(x0)是所有已测站点的第一地区综合系数CP(xi)的线性加权和，即：
Figure BDA0002663415900000101
式中，CP(x0)为待测站点的第二地区综合系数，CP(xi)为第i个已测站点的第一地区综合系数，
Figure BDA0002663415900000102
为待测站点的属性坐标系，
Figure BDA0002663415900000103
为第i个已测站点的属性坐标系，B1,B2,…,Bk为普通克里金法因子，λi为第i个已测站点的克里金权重，n为已测站点的个数，k为普通克里金法因子个数。
进一步的，所述第i个已测站点的克里金权重λi的计算公式为：
Figure BDA0002663415900000104
式中，
Figure BDA0002663415900000105
为第i个和第j个已测站点xi与xj之间的半方差，
Figure BDA0002663415900000106
为第j个已测站点xi与待测站点x0之间的半方差，μ为拉格朗日乘子。
(五)水文设计值预测
步骤5、根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，并基于所述区域回归方程计算待测站点指定设计频率下的第二水文设计值。
具体而言，根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，代入到步骤3建立的区域回归方程中，即可计算得到待测站点指定设计频率下的第二水文设计值，该第二水文设计值为初始值。对于已测站点和待测站点来说，所述区域回归方程的回归方程形式和回归方程参数均相等。
(六)水文设计值修正
在本实施例中，还可以包括：
步骤6、建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数，分别计算各已测站点的第一综合系数的预测值；
步骤7、根据各已测站点的第一地区综合系数的预测值、对应的回归法因子以及回归方程参数，分别计算各已测站点的第一水文设计值的预测值；
步骤8、根据各已测站点的第一水文设计值及其预测值计算第一水文设计值的误差值，建立基于普通克里金法因子的反属性插值方法，根据所有已测站点的第一水文设计值的误差值计算待测站点的第二水文设计值的误差值；
具体而言，根据步骤2中计算的已测站点指定设计频率下的第一水文设计值，以及步骤7计算的已测站点指定频率下的第一水文设计值的预测值，计算两者的误差值，其计算公式为：
Figure BDA0002663415900000111
式中，
Figure BDA0002663415900000112
为第i个已测站点指定设计频率下的第一水文设计值，
Figure BDA0002663415900000113
为第i个已测站点指定设计频率下的第一水文设计值的预测值，i为第i个已测站点，n为已测站点的个数。
所述建立基于普通克里金法因子的反属性插值方法，根据所有已测站点的第一水文设计值的误差值计算待测站点的第二水文设计值的误差值
Figure BDA0002663415900000114
其计算公式为：
Figure BDA0002663415900000115
其中，
Figure BDA0002663415900000116
式中，
Figure BDA0002663415900000117
为第i个已测站点指定设计频率下的第一水文设计值的误差值，di为待测站点与第i个已测站点之间的属性距离，
Figure BDA0002663415900000118
为待测站点的属性坐标系，
Figure BDA0002663415900000119
Figure BDA00026634159000001110
为第i个已测站点的属性坐标系，B1,B2,…,Bk为普通克里金法因子，n为已测站点的个数，k为普通克里金法因子个数。
步骤9、根据所述第二水文设计值的误差值计算待测站点的第二水文设计值的修正值。
将步骤5计算得到的待测站点指定设计频率下的第二水文设计值和步骤8计算得到的第二水文设计值的误差值相加，即可得到待测站点指定设计频率下的第二水文设计值的修正值，进而可明显提高水文设计值计算结果的准确性。
(七)水文设计值优化
步骤10、重复步骤1至9，计算待测站点所有设计频率下的第二水文设计值的修正值；根据水利水电工程水文设计实际情况，设计频率P常取：
P＝0.1％、0.2％、0.5％、1％、2％、5％、10％、20％、50％、75％、90％、95％。
步骤11、建立基于SCE-UA算法的水文频率分析方法，计算待测站点的水文频率分布参数；
假定水文变量X服从P-III型分布，将待测站点所有设计频率下的第二水文设计值的修正值及其对应的设计频率作为样本点据；
建立基于SCE-UA算法的水文频率分析方法，计算水文频率分布参数，分别为期望值Ex、变差系数Cv、偏态系数Cs。
步骤12、由水文频率分布线型及其参数值，计算待测站点所有设计频率下的第二水文设计值的最终值。
具体的，根据水文频率分布线型及其参数值，计算待测站点所有设计频率下的第二水文设计值，该第二水文设计值为最终值，其计算公式为：
QP＝Ex×(1+ΦP×Cv)；
其中，
Figure BDA0002663415900000121
式中，P为指定设计频率，Ex为期望值，Cv为变差系数，Cs为偏态系数，ΦP为离均系数，Γ为Gamma函数。
(八)方法可靠性评估
步骤13、针对指定设计频率下，从所有已测站点中抽取一站作为验证站，剩余已测站点组合作为率定站；
具体的，采用交叉验证中的留一法，每次从所有已测站点中只抽取1站作为验证站，进行方法预测，并假设该站无水文设计值成果；
剩余已测站点组合作为率定站，进行方法训练。
步骤14、假设率定站和验证站分别为已测站点和待测站点，根据率定站指定设计频率下的第一水文设计值，计算验证站指定设计频率下的第二水文设计值的预测值；
步骤15、重复步骤13至14，将所有已测站点逐一抽取作为验证站，并将验证站假设为待测站点；依次取不同的设计频率，计算所有验证站所有设计频率下的第二水文设计值的预测值，并将其与对应设计频率下的第二水文设计值进行比较，进而评估该方法的可靠性。
具体的，本实施例优先采用相对误差(Relative Error，)、平均相对误差绝对值(Absolute of Mean Relative Error，AMRE)来评估该方法的可靠性；其中，RE用来评估方法的无偏性，其值越小越好，AMRE用来评估方法的有效性，其值越小越好；其计算公式分别为：
Figure BDA0002663415900000122
Figure BDA0002663415900000123
式中，
Figure BDA0002663415900000124
为第i个已测站点指定频率下的水文设计值的预测值，
Figure BDA0002663415900000125
为第i个已测站点指定频率下的水文设计值，i为第i个已测站点，n为已测站点的个数，k为第k个设计频率，p为指定设计频率的个数。

Claims (10)

1.一种基于空间统计理论的水文设计值计算方法，其特征在于，包括以下步骤：
步骤1、基于可获取的多源海量遥感大数据，提取研究流域及其邻近流域的多种属性栅格数据；建立基于地理探测器的属性因子优选方法，从所述属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集；
步骤2、收集并整理多个已测站点的水文资料，根据所述水文资料并建立基于SCE-UA算法的水文频率分析方法，计算各已测站点指定设计频率下的第一水文设计值；
步骤3、对于各已测站点，采用逐步回归法，从所述关键属性因子集中选取对第一水文设计值影响显著的属性因子作为回归法因子；根据所有已测站点的第一水文设计值和对应的回归法因子，建立区域回归方程，计算各已测站点的第一地区综合系数以及回归方程参数；
步骤4、对于各已测站点，采用逐步普通克里金法，从所述关键属性因子集中选取对第一地区综合系数影响显著的属性因子作为普通克里金法因子；建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数；
步骤5、根据待测站点的第二地区综合系数、对应的回归法因子以及回归方程参数，并基于所述区域回归方程计算待测站点指定设计频率下的第二水文设计值。
2.如权利要求1所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，还包括：
步骤6、建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数，分别计算各已测站点的第一综合系数的预测值；
步骤7、根据各已测站点的第一地区综合系数的预测值、对应的回归法因子以及回归方程参数，分别计算各已测站点的第一水文设计值的预测值；
步骤8、根据各已测站点的第一水文设计值及其预测值计算第一水文设计值的误差值，建立基于普通克里金法因子的反属性插值方法，根据所有已测站点的第一水文设计值的误差值计算待测站点的第二水文设计值的误差值；
步骤9、根据所述第二水文设计值的误差值计算待测站点的第二水文设计值的修正值。
3.如权利要求2所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，还包括：
步骤10、重复步骤1至9，计算待测站点所有设计频率下的第二水文设计值的修正值；
步骤11、建立基于SCE-UA算法的水文频率分析方法，计算待测站点的水文频率分布参数；
步骤12、由水文频率分布线型及其参数值，计算待测站点所有设计频率下的第二水文设计值的最终值。
4.如权利要求2所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，还包括：
步骤13、针对指定设计频率下，从所有已测站点中抽取一站作为验证站，剩余已测站点组合作为率定站；
步骤14、假设率定站和验证站分别为已测站点和待测站点，根据率定站指定设计频率下的第一水文设计值，计算验证站指定设计频率下的第二水文设计值的预测值；
步骤15、重复步骤13至14，将所有已测站点逐一抽取作为验证站，并将验证站假设为待测站点；依次取不同的设计频率，计算所有验证站所有设计频率下的第二水文设计值的预测值，并将其与对应设计频率下的第二水文设计值进行比较，进而评估该方法的可靠性。
5.如权利要求2所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，步骤8中，所述第二水文设计值的误差值
Figure FDA0002663415890000021
的计算公式为：
Figure FDA0002663415890000022
其中，
Figure FDA0002663415890000023
式中，
Figure FDA0002663415890000024
为第i个已测站点指定设计频率下的第一水文设计值的误差值，di为待测站点与第i个已测站点之间的属性距离，
Figure FDA0002663415890000025
为待测站点的属性坐标系，
Figure FDA0002663415890000026
Figure FDA0002663415890000027
为第i个已测站点的属性坐标系，B1，B2，…，Bk为普通克里金法因子，n为已测站点的个数，k为普通克里金法因子个数。
6.如权利要求1所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，步骤1中，所述从属性栅格数据中选取决定水文变量空间分布规律的关键属性因子集具体包括：
采用地理探测器中的因子探测法，将水文变量栅格数据作为因变量，其余属性栅格数据作为自变量，探测因变量的空间分异性以及各自变量多大程度上解释了因变量的空间分异，选取对水文变量影响显著的属性组成关键属性因子集。
7.如权利要求1所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，步骤2中，所述计算各已测站点指定设计频率下的第一水文设计值具体包括：
假定水文变量X服从P-III型分布，记作X～Γ(x；a，α，β)，其概率密度函数f(x)和分布函数F(x)计算公式为：
Figure FDA0002663415890000028
式中，a为位置参数，α为形状参数，β为尺度参数；
根据指定的水文频率分布线型即P-III型分布，建立基于SCE-UA算法的水文频率分析方法，计算水文频率分布参数；
根据水文频率分布线型及其参数值，计算各已测站点指定设计频率下的第一水文设计值QP，其计算公式为：
QP＝Ex×(1+ΦP×Cv)；
其中，
Figure FDA0002663415890000031
式中，P为设计频率，Ex为期望值，Cv为变差系数，Cs为偏态系数，ΦP为离均系数，Γ为Gamma函数。
8.如权利要求1所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，步骤3中，所述区域回归方程采用多元幂函数方程，其计算公式为：
Figure FDA0002663415890000032
式中，QP为指定设计频率下的水文设计值，CP为地区综合系数，A1，A2，…，Am为回归法因子，α1，α2，…，αm为回归方程参数，m为回归法因子个数。
9.如权利要求1所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，步骤4中，所述建立基于属性因子的普通克里金法，根据所有已测站点的第一地区综合系数计算待测站点的第二地区综合系数具体包括：
假设各已测站点的第一地区综合系数分别为CP(xi)(i＝1，2，…，n)，则待测站点的第二地区综合系数CP(x0)是所有已测站点的第一地区综合系数CP(xi)的线性加权和，即：
Figure FDA0002663415890000033
式中，CP(x0)为待测站点的第二地区综合系数，CP(xi)为第i个已测站点的第一地区综合系数，
Figure FDA0002663415890000034
为待测站点的属性坐标系，
Figure FDA0002663415890000035
为第i个已测站点的属性坐标系，B1，B2，…，Bk为普通克里金法因子，λi为第i个已测站点的克里金权重，n为已测站点的个数，k为普通克里金法因子个数。
10.如权利要求9所述的一种基于空间统计理论的水文设计值计算方法，其特征在于，所述第i个已测站点的克里金权重λi的计算公式为：
Figure FDA0002663415890000036
式中，
Figure FDA0002663415890000037
为第i个和第j个已测站点xi与xj之间的半方差，
Figure FDA0002663415890000038
为第j个已测站点xi与待测站点x0之间的半方差，μ为拉格朗日乘子。
CN202010911383.3A 2020-09-02 2020-09-02 一种基于空间统计理论的水文设计值计算方法 Active CN112035553B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010911383.3A CN112035553B (zh) 2020-09-02 2020-09-02 一种基于空间统计理论的水文设计值计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010911383.3A CN112035553B (zh) 2020-09-02 2020-09-02 一种基于空间统计理论的水文设计值计算方法

Publications (2)

Publication Number Publication Date
CN112035553A CN112035553A (zh) 2020-12-04
CN112035553B true CN112035553B (zh) 2022-11-29

Family

ID=73591286

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010911383.3A Active CN112035553B (zh) 2020-09-02 2020-09-02 一种基于空间统计理论的水文设计值计算方法

Country Status (1)

Country Link
CN (1) CN112035553B (zh)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109523066A (zh) * 2018-10-29 2019-03-26 东华理工大学 一种基于克里金插值的pm2.5新增移动站点选址方法
CN111414671A (zh) * 2019-12-20 2020-07-14 福建省水利水电勘测设计研究院 小流域暴雨推求洪水的方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109523066A (zh) * 2018-10-29 2019-03-26 东华理工大学 一种基于克里金插值的pm2.5新增移动站点选址方法
CN111414671A (zh) * 2019-12-20 2020-07-14 福建省水利水电勘测设计研究院 小流域暴雨推求洪水的方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于DEM的榆林市降水空间插值方法分析;刘智勇等;《西北农林科技大学学报(自然科学版)》;20100715(第07期);全文 *

Also Published As

Publication number Publication date
CN112035553A (zh) 2020-12-04

Similar Documents

Publication Publication Date Title
Wang et al. Predictive mapping of soil total nitrogen at a regional scale: A comparison between geographically weighted regression and cokriging
Masih et al. Assessing the impact of areal precipitation input on streamflow simulations using the SWAT Model 1
Vogel Regional calibration of watershed models
Wang et al. Comparison of geographically weighted regression and regression kriging for estimating the spatial distribution of soil organic matter
CN103093114A (zh) 一种基于地形和土壤特性的分布式流域缺水量测算方法
CN111985389B (zh) 一种基于流域属性距离的流域相似判别方法
Chea et al. Flow simulation in an ungauged catchment of Tonle Sap Lake Basin in Cambodia: Application of the HEC-HMS model
CN105528523A (zh) 一种基于遥感数据的土壤厚度反演方法
CN117408495A (zh) 一种基于土地资源综合管理的数据分析方法及系统
CN113591288A (zh) 一种基于克里格插值的土壤湿度数据预测方法及装置
CN114461983B (zh) 一种基于水量平衡原理的卫星降水产品空间降尺度方法
Wolff et al. Toward geostatistical unbiased predictions of flow duration curves at ungauged basins
Zhang et al. Spatial patterns and controlling factors of the evolution process of karst depressions in Guizhou province, China
Nourani et al. Spatiotemporal assessment of groundwater quality and quantity using geostatistical and ensemble artificial intelligence tools
CN113762615A (zh) 洪水预测方法、装置、计算机设备和存储介质
CN112035553B (zh) 一种基于空间统计理论的水文设计值计算方法
CN117494042A (zh) 一种基于气候相似性迁移的缺资料地区土壤水融合方法
Wijemannage et al. Comparison of spatial interpolation methods for rainfall data over Sri Lanka
Tu et al. Spatial downscaling analysis of GPM IMERG precipitation dataset based on multiscale geographically weighted regression model: a case study of the Inner Mongolia Reach of the Yellow River basin
Gadgay et al. Novel ensemble neural network models for better prediction using variable input approach
CN114385959A (zh) 一种近坝区子流域单元划分方法、装置及存储介质
Feng et al. Predicting soil depth in a large and complex area using machine learning and environmental correlations
Saha et al. An integrated multi-criteria decision analysis and geographic information system-based assessment of groundwater potentiality and stress zones for sustainable agricultural practices: a case study of agriculture-dominating Koch Bihar District, West Bengal
İçağa et al. Comparative analysis of different interpolation methods in modeling spatial distribution of monthly precipitation
Li et al. Improvement of the multi-source weighted-ensemble precipitation dataset and application in the arid area of Tianshan Mountains, central Asia

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant