CN115796054B - 一种情景发现和脆弱性分析方法、系统、终端及存储介质 - Google Patents

Abstract

本发明涉及情景发现和脆弱性分析方法、系统、终端及存储介质，方法包括以下步骤：在采用其他外部仿真模型得到对情景集的模拟产生的案例数据库后，判断案例的数量是否超过阈值，超出则采用PRIM算法直接计算获得关键政策情景与脆弱性分析结果，不超出则采用帕累托最优集算法进行求解获得最优集，如果最优集中案例数量超过阈值，则再对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，如果最优集中案例数量不超过阈值，则将最优集作为关键政策情景与脆弱性分析结果；采用本发明的结合方式能高效且简明的完成稳健决策中“情景发现和脆弱性分析”这一步骤。

Description

一种情景发现和脆弱性分析方法、系统、终端及存储介质

技术领域

本发明涉及水资源等气候相关技术领域，更具体地说，涉及一种情景发现和脆弱性分析方法、系统、终端及存储介质。

背景技术

在传统的决策过程中，对政策的脆弱性分析通常采用情景分析法。传统的情景分析法在面对简单而确定的未来时有较好的表现。但当决策涉及到多种措施、不确定因素和指标时，对情景集的模拟将产生成千上万的案例。在这种情况下，传统的情景分析法无法进行定量决策分析，同时又因为情景案例数量过多而难以实施。此时，先进的数据分析技术可以提供帮助。“情景发现”即使用统计或数据挖掘算法，在海量的多因子仿真模拟模型结果的数据库中找到与决策政策相关的案例集，发现关键情景并转化为易于解释的不确定因子的组合，且这些因子组合对政策相关的情况有很强的预测能力。找出的案例集能够被方便地解释为不同的未来情景，有助于阐明和量化在深度不确定情况下各未来决策路径之间的比较与权衡。此外，情景发现还可以帮助用户确定他们所将要面临的最重要的情景。

目前较为常用的工程技术方法包括Patient Rule Induction Method（PRIM）算法和采用各类算法求帕累托最优集等。由美国RAND公司开发的PRIM算法是基于R语言的，无法编译为可独立运行的可执行文件，且并不提供可互动的可视化界面，需要使用者自行编写代码实现可视化，操作繁琐，对于非专业使用者并不友好。现有的求多目标优化的帕累托最优集的算法中，较为先进的是非支配排序遗传算法（NSGA）。在传统设计行业中，NSGA算法通常被设计用于解决可用数学公式描述的多目标决策问题，而在雨水利用分析和洪涝风险识别等水相关领域中，具体的决策问题通常无法使用简单的数学公式进行描述。一般来说，这些水相关问题都是由复杂的水文-水动力-水质关系模型进行模拟，从模拟结果中由决策者选取出特定的多个效益指标进行综合评估。此时，NSGA算法需要解决的问题是从多维空间的离散数据点中求出离散的帕累托最优集，这与连续的理论问题有较大区别。因此，常见的用于理论问题的NSGA算法并不能完全满足合雨水利用分析和洪涝风险识别决策的相关需求。

综上所述，对于雨水利用分析和洪涝风险识别等水资源利用和规划的相关领域的稳健决策，决策者若要使用现有的常见算法，需要拥有丰富的专业知识和足够的实践经验以解决决策实施环节中可能产生的问题。在所有的决策环节中，情景发现和脆弱性分析对决策结果至关重要，是决策者通常需要特别关心和仔细进行分析的环节。对非专业的决策者而言，直接采用现有的常见算法进行这一环节的分析在理论和实际应用上的学习成本和试错成本较高。因此，需要一种能够解决该问题的更为友好并易于理解和使用的情景发现和脆弱性分析方式方法。

发明内容

本发明要解决的技术问题在于，针对现有技术的上述缺陷，提供一种情景发现和脆弱性分析方法，还提供了一种情景发现和脆弱性分析系统、一种情景发现和脆弱性分析终端及一种计算机可读存储介质。

本发明解决其技术问题所采用的技术方案是：

构造一种情景发现和脆弱性分析方法，其中，包括以下步骤：

判断案例数据库中的对情景集的模拟产生的案例的数量是否超过第一阈值；

是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；

否则采用帕累托最优集算法进行求解获得最优集后，判断最优集中案例数量是否超过第二阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果。

本发明所述的情景发现和脆弱性分析方法，其中，案例的参数包括案例对应的情景编号ID、情景设置中各因子变量和案例仿真模拟的多维目标的绩效结果。

本发明所述的情景发现和脆弱性分析方法，其中，所述PRIM算法实现方式包括：

第一步：初始化参数，参数包括目标绩效阈值与数据绩效和阈值关系的大小关系判断方法、控制目标函数、剔除参数alpha_1、粘贴参数alpha_2和最小质量mass_min；将决策问题抽象为f(y|x)，根据因子变量{X1, X2, …, Xn}，预测满足多维目标{Y1, Y2, …,Yn}中某一个目标发生的概率P(Y=1|X)；

第二步：求所有因子变量当前在alpha_1和1-alpha_1分位数的值，即统计当前数据点在各因子变量的累积分布概率等于alpha_1和1-alpha_1处的值，得到X（alpha_1）和X（1-alpha_1）；

第三步：选取各因子变量在分位数0至alpha_1、1-alpha_1至1之间的区域，生成备选的剔除空间区域；

第四步：判断是否有任一备选区域的剔除能使进行剔除后得到的空间区域的控制目标函数的结果比未进行剔除的原始空间区域的控制目标函数的结果大，是则进行下一步，否则跳转至第七步；

第五步：剔除能让目标函数结果增大最多的备选区域；

第六步：判断剩余区域形成新的超立方体的质量是否大于最小质量mass_min，是则进行下一步，否则跳转至第二步；

第七步：求所有因子变量当前向已剔除的区域往外扩展alpha_2分位数的值，即0-alpha_2和1+alpha_2分位数处的值，得到X（0-alpha_2）和X（1+alpha_2）；

第八步：选取各因子变量在分位数0-alpha_2至0、1至1+alpha_2之间的区域，生成备选的粘贴空间区域；

第九步：判断是否有任一备选区域的粘贴能使进行粘贴后得到的空间区域的控制目标函数的结果比未进行粘贴的原始空间区域的控制目标函数的结果大，是则粘贴让目标函数结果增大最多的备选区域并跳转至第七步，否则进行下一步；

第十步：输出案例数据库的目标规则集和因子变量关键变化范围；目标规则集包含控制目标函数值随各因子变量区域的剔除和粘贴的变化规律与特点，因子变量关键变化范围为最后剩余区域的变量上下限区间。

本发明所述的情景发现和脆弱性分析方法，其中，所述第二步至第六步采用方法：

控制目标函数为：；

其中，表示在m维空间的当前超立方体里的样本量，为案例数据库中的绩效结果，f(y)为当前超立方体里的正样本浓度，当前超立方体在各个维度上的上下限与样本群在该维度上的自变量上下限相同；

剔除参数alpha_1表示每次剔除的样本比例，低于分位数X（alpha_1）或高于分位数X（1-alpha_1）的样本将会被剔除；

依据剔除参数alpha_1，在当前超立方体内部生成一组备选的剔除空间区域：

其中：

；

；

和分别表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体里的α和1-α分位数的值；从所有的备选区域中选择一个，使得剔除该区域后，剩余空间区域的控制目标函数的结果比剔除前增加最多；从剔除一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选剔除区域为：

；

不断迭代执行剔除操作，直到备选区域集合中所有备选区域的剔除都不能增大PRIM控制目标函数结果，或满足以下条件：

；

指示函数的含义为：如果样本在当前超立方体里，则取值为1，反之为0，n是总样本量；是属于当前超立方体的样本量相对于总样本量的占比，即PRIM过程中超立方体的质量；是目标超立方体的最小样本占比，即预设的PRIM过程中超立方体的最小质量mass_min；当超立方体的质量达到了用户预先设定的最小质量时，PRIM的迭代剔除步骤结束。

本发明所述的情景发现和脆弱性分析方法，其中，所述第七步至第十步采用方法：

根据粘贴参数alpha_2，在当前超立方体外部生成一组备选区域：

；

其中，

；

；

和表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体向外扩展α分位数的值，和表示Xj在中的最小值和最大值；从所有的备选区域中选择一个，使得粘贴该区域后，剩下区域的控制目标函数的结果比粘贴前增加最多；为粘贴一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选粘贴区域为：

；

不断迭代执行粘贴操作，直到备选区域集合中所有备选区域的粘贴都不能增大PRIM控制目标函数结果；迭代粘贴操作停止后，PRIM输出最终的超立方体作为结果，结果包含了案例数据库的目标规则集和因子变量的关键变化范围信息。

本发明所述的情景发现和脆弱性分析方法，其中，所述帕累托最优集算法包括方法：

对多个目标采用以下公式进行优化：

；

；

；

；

其中o为多维目标F(x)={Y1, Y2, …, Yn}的目标个数，u为不等式约束的个数，v为等式约束的个数，和分别为中第i个变量x_i的上界和下界，w为变量x_i的总个数；

采用遗传算法，求帕累托最优集。

本发明所述的情景发现和脆弱性分析方法，其中，所述遗传算法包括内容：

随机产生初始种群，采用非支配排序通过遗传算法的选择、交叉、变异三个基本操作得到第一代子代种群；

将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离迭代进行拥挤度计算：

；

为个体在目标空间中的曼哈顿距离；

对局部的支配关系排序为：

；

其中，和分别为局部空间中i和j两个个体的排序，和分别为i和j由曼哈顿距离所代表的个体拥挤度；

使用二元锦标赛选择方法或者参考点选择方法来维持种群的多样性；

根据非支配关系以及个体的拥挤度选取合适的个体组成新的父代种群；对于两个样本点，优先选择帕累托前沿排序更小的；如果两个样本点在同一个帕累托前沿上，则选择拥挤度更大的样本点；

判断是否满足预设的结束条件，是则结束，否则调整至将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离进行拥挤度计算步骤。

一种情景发现和脆弱性分析系统，应用于如上述的情景发现和脆弱性分析方法，其中，包括案例数据库、数据处理模块和数据输出模块；

所述案例数据库，存储对情景集的模拟产生的案例；

所述数据处理模块，判断案例数据库中的对情景集的模拟产生的案例的数量是否超过阈值；是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；否则采用帕累托最优集算法进行求解获得最优集后，判断最优集中案例数量是否超过阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果；

所述数据输出模块，将关键政策情景与脆弱性分析结果输出。

一种情景发现和脆弱性分析终端，包括存储器、处理器以及存储在所述存储器中并可在所述处理器上运行的计算机程序，其中，所述处理器执行所述计算机程序时实现如上述方法的步骤。

一种计算机可读存储介质，所述计算机可读存储介质存储有计算机程序，其中，所述计算机程序被处理器执行时实现如上述方法的步骤。

本发明的有益效果在于：在采用其他外部仿真模型得到对情景集的模拟产生的案例数据库后，判断案例的数量是否超过阈值，超出则采用PRIM算法直接计算获得关键政策情景与脆弱性分析结果，不超出则采用帕累托最优集算法进行求解获得最优集，如果最优集中案例数量是否超过阈值，则再对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，如果最优集中案例数量不超过阈值，则将最优集作为关键政策情景与脆弱性分析结果；采用本发明的结合方式能高效且简明的完成稳健决策中“情景发现和脆弱性分析”这一步骤。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将结合附图及实施例对本发明作进一步说明，下面描述中的附图仅仅是本发明的部分实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他附图：

图1是本发明较佳实施例的情景发现和脆弱性分析方法逻辑流程图；

图2是本发明较佳实施例的情景发现和脆弱性分析方法PRIM算法逻辑流程图；

图3是本发明较佳实施例的情景发现和脆弱性分析方法PRIM算法权衡过程示意图；

图4是本发明较佳实施例的情景发现和脆弱性分析方法PRIM过程中超立方体细节示意图；

图5是本发明较佳实施例的情景发现和脆弱性分析方法NSGA2/NGSA3算法流程图；

图6是本发明较佳实施例的情景发现和脆弱性分析方法NGSA2的拥挤度计算示意图；

图7是本发明较佳实施例的情景发现和脆弱性分析方法NGSA2的二元锦标赛选择方法示意图；

图8是本发明较佳实施例的情景发现和脆弱性分析方法NSGA3的前沿划分示意图；

图9是本发明较佳实施例的情景发现和脆弱性分析方法NGSA3的参考点选择方法示意图；

图10是本发明较佳实施例的情景发现和脆弱性分析方法一个多维帕累托优化集在Y1-Y2平面上的投影示例；

图11是本发明较佳实施例的情景发现和脆弱性分析系统原理框图。

具体实施方式

为了使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的技术方案进行清楚、完整的描述，显然，所描述的实施例是本发明的部分实施例，而不是全部实施例。基于本发明的实施例，本领域普通技术人员在没有付出创造性劳动的前提下所获得的所有其他实施例，都属于本发明的保护范围。

本发明较佳实施例的情景发现和脆弱性分析方法，如图1所示，同时参阅图2-10，方法包括以下步骤：

判断案例数据库中的对情景集的模拟产生的案例的数量是否超过第一阈值；是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；否则采用帕累托最优集算法进行求解获得最优集后，判断最优集中案例数量是否超过第二阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果；

在实际的工程决策中，基于案例数据库，决策者需要进行关键情景发现和脆弱性分析，需要考虑在怎样的不确定条件下，所有决策都不能满足利益相关者的目标？在什么样的未来情景中，一个特定的措施组合会失败？为了辅助非专业使用者进行稳健决策中的“情景发现和脆弱性分析”这一环节，使有经验的用户在更方便、更准确地实施稳健决策框架，本专利集成了PRIM方法和基于离散数据求多目标帕累托最优集两个算法模块（如图1所示）；

在采用其他外部仿真模型得到模拟结果的数据库后，用户应根据自身分析需求和实际研究情况，选择其中一种算法。这里所说的结果数据库由不同情景的仿真模拟结果所组成，每一个案例数据包含三部分：案例对应的情景编号ID、情景设置中各因子变量{X1,X2, …, Xn}和案例仿真模拟的多维目标的绩效结果{Y1, Y2, …, Yn}。一般来说，当数据库中的案例数量较少时，可以采用基于离散数据求多目标帕累托最优集的算法，直接得到关键政策情景与其脆弱性；当数据库中的案例数量较多（上千）时，用户应采用PRIM方法确保情景发现和脆弱性分析的过程简洁明了且易于解释。用户也可以先求出数据库的帕累托最优集，若该最优集内关键政策情景数量太多，再采用PRIM进行分析。两种方法得到结果后，用户需结合具体的研究问题和利益相关方的需求等实际情况，完成稳健决策中“情景发现和脆弱性分析”这一步骤。

1. PRIM方法（如图2所示）

当决策涉及到多种措施、不确定因素和指标时，对情景集的模拟将产生成千上万的案例。当有两个以上的不确定因子时，很难通过简单的方法来识别脆弱情景。PRIM方法有助于找到并简要地解释措施失败或成功的未来条件，其基本理念是将稳健决策框架中对未来情景模拟结果的评估产生的大量数据提炼为用户和利益相关者最关心的信息。当有大量的未来情景案例需要分析时，PRIM可用于确定哪些不确定因子在系统易受影响的未来情景中最为重要，以及定义一套简明的规则以描述不确定因子的变化范围。

将决策问题抽象为f(y|x)，根据因子变量{X1, X2, …, Xn}，预测满足多维目标{Y1, Y2, …, Yn}中某一个目标发生的概率P(Y=1|X)。但我们其实更关心满足条件的Y值案例的浓度很高（或很低）的某个局部空间。现在若有N个样本{Xi,Y}^N，我们希望从M维变量空间寻找一个子空间，使得这个子空间的目标变量浓度尽可能高。这个问题被称为子群识别。常见的操作是在M维空间里横竖切几刀，剔除不重要的区域，最后找到并保留目标区域。简单来说，PRIM是一种多维度的规则发现方法，基本的操作过程就是自上向下削除样本空间。

PRIM算法实现方式包括：

第一步：初始化参数，参数包括目标绩效阈值与数据绩效和阈值关系的大小关系判断方法、控制目标函数、剔除参数alpha_1、粘贴参数alpha_2和最小质量mass_min；将决策问题抽象为f(y|x)，根据因子变量{X1, X2, …, Xn}，预测满足多维目标{Y1, Y2, …,Yn}中某一个目标发生的概率P(Y=1|X)；

第二步：求所有因子变量当前在alpha_1和1-alpha_1分位数的值，即统计当前数据点在各因子变量的累积分布概率等于alpha_1和1-alpha_1处的值，得到X（alpha_1）和X（1-alpha_1）；

第三步：选取各因子变量在分位数0至alpha_1、1-alpha_1至1之间的区域，生成备选的剔除空间区域；

第四步：判断是否有任一备选区域的剔除能使进行剔除后得到的空间区域的控制目标函数的结果比未进行剔除的原始空间区域的控制目标函数的结果大，是则进行下一步，否则跳转至第七步；

第五步：剔除能让目标函数结果增大最多的备选区域；

第六步：判断剩余区域形成新的超立方体的质量是否大于最小质量mass_min，是则进行下一步，否则跳转至第二步；

第七步：求所有因子变量当前向已剔除的区域往外扩展alpha_2分位数的值，即0-alpha_2和1+alpha_2分位数处的值，得到X（0-alpha_2）和X（1+alpha_2）；

第八步：选取各因子变量在分位数0-alpha_2至0、1至1+alpha_2之间的区域，生成备选的粘贴空间区域；

第九步：判断是否有任一备选区域的粘贴能使进行粘贴后得到的空间区域的控制目标函数的结果比未进行粘贴的原始空间区域的控制目标函数的结果大，是则粘贴让目标函数结果增大最多的备选区域并跳转至第七步，否则进行下一步；

第十步：输出案例数据库的目标规则集和因子变量关键变化范围；目标规则集包含控制目标函数值随各因子变量区域的剔除和粘贴的变化规律与特点，因子变量关键变化范围为最后剩余区域的变量上下限区间。

优选的，第二步至第六步采用方法：

控制目标函数为：

；

其中，表示在m维空间的当前超立方体里的样本量，为案例数据库中的绩效结果，f(y)为当前超立方体里的正样本浓度，当前超立方体在各个维度上的上下限与样本群在该维度上的自变量上下限相同；

剔除参数alpha_1表示每次剔除的样本比例，低于分位数X（alpha_1）或高于分位数X（1-alpha_1）的样本将会被剔除；

依据剔除参数alpha_1，在当前超立方体内部生成一组备选的剔除空间区域：

其中：

；

；

和分别表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体里的α和1-α分位数的值；从所有的备选区域中选择一个，使得剔除该区域后，剩余空间区域的控制目标函数的结果比剔除前增加最多；从剔除一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选剔除区域为：

；

不断迭代执行剔除操作，直到备选区域集合中所有备选区域的剔除都不能增大PRIM控制目标函数结果，或满足以下条件：

；

指示函数的含义为：如果样本在当前超立方体里，则取值为1，反之为0，n是总样本量；是属于当前超立方体的样本量相对于总样本量的占比，即PRIM过程中超立方体的质量；是目标超立方体的最小样本占比，即预设的PRIM过程中超立方体的最小质量mass_min；当超立方体的质量达到了用户预先设定的最小质量时，PRIM的迭代剔除步骤结束。

在剔除步骤中，决策边界是局部最优，依赖于上一步结果。也就是说，剔除是贪婪的。贪婪算法在搜索过程中是局部最优，而非全局最优。粘贴是一种自下向上的操作，是与剔除相反的操作。因此将会上升，而同时目标函数也可能会增大，这一操作试图让PRIM跳出局部而考虑全局最优。

粘贴与剔除操作的流程类似，粘贴的步骤第七步至第十步采用方法：

根据粘贴参数alpha_2，在当前超立方体外部生成一组备选区域：

；

其中，

；

；

和表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体向外扩展α分位数的值，，和表示Xj在中的最小值和最大值；与剔除参数类似，粘贴参数也选用一个比较小的数值较为合适（0.05-0.1）；从所有的备选区域中选择一个，使得粘贴该区域后，剩下区域的控制目标函数的结果比粘贴前增加最多；为粘贴一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选粘贴区域为：

；

不断迭代执行粘贴操作，直到备选区域集合中所有备选区域的粘贴都不能增大PRIM控制目标函数结果；因为粘贴操作不会使新的超立方体的质量减小，所以与剔除操作不同，这里不需要考虑最小质量。迭代粘贴操作停止后，PRIM输出最终的超立方体作为结果，结果包含了案例数据库的目标规则集和因子变量的关键变化范围信息，供用户进行进一步的分析。

综上所述，PRIM算法流程可以总结为：先迭代执行剔除操作，直到所有备选区域的剔除都不能改善超立方体的PRIM控制目标函数，或超立方体的样本数相对于总样本数的占比小于预设的最小质量；再尝试迭代进行粘贴操作，直到所有备选区域的粘贴都不能改善超立方体的PRIM控制目标函数；最后，输出最终的超立方体作为计算结果。

在情景发现和脆弱性分析中，PRIM迭代地确定了一个超立方体，针对多维目标{Y1, Y2, …, Yn}中用户选择的一个目标Y，由包含的不同变量{X1, X2, …, Xn}的范围来定义。在实际使用PRIM时，我们需要知道如何确定PRIM何时找到了正确描述脆弱性的条件。一个特定的脆弱性条件是否有效，在PRIM中可以用以下三个特征来衡量：

（1）覆盖率：脆弱条件捕捉到的决策失败案例占所有失败案例的比例是多少？理想情况下，脆弱条件可以描述所有的决策失败案例，即覆盖率为100%。

（2）密度：在脆弱条件描述的所有案例中，决策失败案例所占的比例是多少？理想情况下，脆弱性中的所有案例都是决策失败的，即密度为100%。

（3）可解释性：用户理解脆弱情景的难度如何？用于定义情景的不确定因子的数量可在某种程度上代表可解释性。一般来说，使用的不确定因子数量越少，可解释性就越高。

我们希望可以平衡覆盖率、密度和可解释性这三个样本集属性。为此，PRIM先迭代剔除区域，缩小各不确定变量的范围，再尝试粘贴，增大各不确定变量范围，目标就是获得更好的在覆盖率和密度之间的权衡。

本发明将PRIM所有超立方体的情况与每个超立方体的细节进行了可视化，以便用户在PRIM迭代中权衡覆盖率和密度，更直观地观察PRIM的过程并选择符合期望的超立方体。每个点代表了一个超立方体，从左上角到右下角为PRIM超立方体逐渐减小和/或增大的过程（图3）。其中横坐标轴为超立方体的密度，纵坐标为超立方体的覆盖率，点的颜色表示所有不确定变量中与超立方体大小相关的变量数目。可以看到，随着PRIM的进行，与超立方体大小相关的不确定变量数目逐渐增多，覆盖率逐渐减小，密度逐渐增大。用户可以将鼠标移动到点上，将会显示超立方体的序号（覆盖率和质量都为100%的初始超立方体为#0）、覆盖率、密度与质量。质量代表了超立方体里的案例数量占空间内案例总数量的比例。

用户还可以用左键点击的特定的超立方体的代表点，对应的超立方体的细节将会显示，便于用户更好的确定每个超立方体的细节（图4）。超立方体的细节图包括了超立方体的具体情况（表现为超立方体超立方体在各个平面上的投影）、超立方体的各个变量的目标区间、表示了该超立方体的统计属性（覆盖率、密度、质量、超立方体内点的结果值平均值、受限制的变量数量）。分析者还可以用右下角的按钮直接转到上一个或下一个超立方体。

为了进行PRIM，用户除了输入各不确定变量{X1, X2, …, Xn}和得到的多维目标{Y1, Y2, …, Yn}中选择的目标Y的值以外，还需要确定并输入PRIM过程的相关参数。这些参数包括：感兴趣案例的阈值、确定感兴趣案例的方法、PRIM评价超立方体质量的目标函数、控制剔除/粘贴过程的参数、超立方体的最小样本占比。

在执行PRIM后，用户一般应该选择最后一个超立方体作为输出结果，但也可以通过整个PRIM权衡过程和各个超立方体的细节，从而主观确定一个超立方体作为结果。PRIM的输出包括了该超立方体的各变量{X1, X2, …, Xn}是否关键和各变量{X1, X2, …,Xn}对特定目标Y的关键变化区间两种信息。基于这些信息，用户应尝试对利益相关方解释不同变量对未来情景和政策措施的影响，以及这些变量在何时会导致政策变得脆弱，而并不专注于数据库中的某一个案例点。用户也可以对多维目标{Y1, Y2, …, Yn}中的每一个目标Yi都进行一次PRIM操作来分析变量{X1, X2, …, Xn}对不同目标的关键变化区间之间的异同。更进一步，用户需要思考哪些特定变量值所组合而成的特定未来情景会是潜在的关键情景，针对这些关键情景，我们应该采取什么政策措施。

基于离散数据求多目标帕累托最优集：

在多目标决策的权衡中，我们需要处理多维目标{Y1, Y2, …, Yn}，即同时优化所有目标Yi。除了PRIM外，另一种较为普适的用于多目标决策的方法是求得多个目标共同的帕累托最优集。对多个目标采用以下公式进行优化：

；

其中o为多维目标F(x)={Y1, Y2, …, Yn}的目标个数，u为不等式约束的个数，v为等式约束的个数，和分别为中第i个变量x_i的上界和下界，w为变量x_i的总个数；；如果一个解决方案在所有目标{Y1, Y2, …, Yn}中具有相等且至少一个更好的结果，那么它就比另一个解决方案更好，我们称为一个解支配了另一个解。如果这两个解互不支配，那么它们就都称为帕累托最优解。对于每个多目标问题，都存在一组帕累托最优解，代表了多个目标{Y1, Y2, …, Yn}之间的最佳权衡，它们的集合称为帕累托最优集。

为了求帕累托最优集，通常采用遗传算法，NSGA2与NSGA3应用较广（如图5所示）。以NSGA2算法为例，NSGA2算法遵循一般的遗传算法的交配和生存选择的步骤：

随机产生初始种群，采用非支配排序通过遗传算法的选择、交叉、变异三个基本操作得到第一代子代种群；

将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离迭代进行拥挤度计算（如图6所示）：

；

为个体在目标空间中的曼哈顿距离；

对局部的支配关系排序为：

；

其中，和分别为局部空间中i和j两个个体的排序，和分别为i和j由曼哈顿距离所代表的个体拥挤度；

为了增加选择压力，使用二元锦标赛选择方法（如图7所示）或者参考点选择方法来维持种群的多样性；

根据非支配关系以及个体的拥挤度选取合适的个体组成新的父代种群；对于两个样本点，优先选择帕累托前沿排序更小的；如果两个样本点在同一个帕累托前沿上，则选择拥挤度更大的样本点；

判断是否满足预设的结束条件，是则结束，否则调整至将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离进行拥挤度计算步骤。

NSGA3与NSGA2具有类似的框架，区别主要在于选择机制的改变，NSGA3对拥挤度排序进行了改进，通过引入参考点来维持种群的多样性，即在初始化算法时在空间上设定一系列参考角，将帕累托前沿进行等距离划分，再按每个划分的前沿进行求解（如图7所示）。如果一个划分的前沿中没有分配任何解，就选择归一化目标空间中垂直距离最小的解。如果添加了此前沿的第二个解决方案，则会随机分配。当算法收敛时，每条前沿都将试图找到一个具有代表性的非支配解。一般来说，NSGA3在优化目标{Y1, Y2, …, Yn}的个数n大于5时对NSGA2有较明显的改善；

在传统的设计行业中，遗传算法主要用来求解多个目标函数，即{Y1, Y2, …,Yn}=f（{X1, X2, …, Xn}），和若干个限制函数所组成的多目标优化问题的帕累托最优集。而在水相关领域中，这些关系模型非常复杂且难以直接通过函数进行描述，需要由复杂的水文-水动力-水质模型进行模拟，再从模拟结果中由用户选取出一些指标Yi作为优化目标。此时需要解决的问题变为从模型模拟结果数据库的离散数据点中求出离散从{X1, X2,…, Xn}映射至{Y1, Y2, …, Yn}的帕累托最优集。

为此，本发明改变了传统的NSGA2和NSGA3的输入模式，能够以数组形式直接读入关系模型模拟得到的{Y1, Y2, …, Yn}指标结果数据集，再将非支配排序的遗传算法应用于指标数据集进行帕累托最优集的求解。与传统的算法求出的帕累托解不同，基于实际数据的帕累托解在目标空间中的分布是不均匀的（受到实际数据点的分布影响）。但在这种情况下，NSGA3的参考点排序会强制使帕累托解的分布变得均匀，在参考点外的帕累托解可能会被剔除（如图8所示）。因此，当本软件中的NSGA2的性能可以接受时，建议优先采用NSGA2求帕累托最优集。

采用NSGA2/NSGA3算法求得数据库中所有案例的帕累托最优集中包括了案例编号ID和案例各目标{Y1, Y2, …, Yn}结果值两种信息。最优集内的每一个案例都代表了一种独特的未来情景和政策措施的结合（如图9所示）。用户应通过案例的编号ID便可以溯源查询对应的情景和政策措施设置，以及因子变量{X1, X2, …, Xn}的取值，并将这些关键案例与利益相关方讨论，从中便可确定一个或多个最佳的政策路径。更进一步，用户可以利用具体的情景的政策组合，来尝试分析解释各不确定性因子{X1, X2, …, Xn}和不同政策措施的脆弱性之间的关联，以及在这些最优集的案例中是否存在某一个或多个由特定的{Y1,Y2, …, Yn}目标取值所决定的关键情景。

本发明可辅助无经验的使用者更快地学习稳健决策，也可协助有相关经验的用户更直观更迅速地获取所需信息，更好地进行框架内各决策步骤，并有可能改善决策分析结果。

一种情景发现和脆弱性分析系统，应用于如上述的情景发现和脆弱性分析方法，如图11所示，包括案例数据库1、数据处理模块2和数据输出模块3；

案例数据库1，存储对情景集的模拟产生的案例；

数据处理模块2，判断案例数据库中的对情景集的模拟产生的案例的数量是否超过阈值；是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；否则采用帕累托最优集算法进行求解获得最优集后，判断最优集中案例数量是否超过阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果；

数据输出模块3，将关键政策情景与脆弱性分析结果输出；

在采用其他外部仿真模型得到对情景集的模拟产生的案例数据库后，判断案例的数量是否超过阈值，超出则采用PRIM算法直接计算获得关键政策情景与脆弱性分析结果，不超出则采用帕累托最优集算法进行求解获得最优集，如果最优集中案例数量是否超过阈值，则再对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，如果最优集中案例数量不超过阈值，则将最优集作为关键政策情景与脆弱性分析结果；采用本发明的结合方式能高效且简明的完成稳健决策中“情景发现和脆弱性分析”这一步骤。

一种情景发现和脆弱性分析终端，包括存储器、处理器以及存储在存储器中并可在处理器上运行的计算机程序，其中，处理器执行计算机程序时实现如上述方法的步骤。

一种计算机可读存储介质，计算机可读存储介质存储有计算机程序，其中，计算机程序被处理器执行时实现如上述方法的步骤。

应当理解的是，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，而所有这些改进和变换都应属于本发明所附权利要求的保护范围。

Claims (10)

1.一种情景发现和脆弱性分析方法，其特征在于，包括以下步骤：

判断案例数据库中的对情景集的模拟产生的案例的数量是否超过第一阈值；

是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；

否则采用帕累托最优集算法对案例数据库进行求解获得最优集后，判断最优集中案例数量是否超过第二阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果；

使用PRIM算法时，将PRIM所有超立方体的情况与每个超立方体的细节进行可视化，供用户在PRIM迭代中权衡覆盖率和密度。

2.根据权利要求1所述的情景发现和脆弱性分析方法，其特征在于，案例的参数包括案例对应的情景编号ID、情景设置中各因子变量和案例仿真模拟的多维目标的绩效结果。

3.根据权利要求2所述的情景发现和脆弱性分析方法，其特征在于，所述PRIM算法实现方式包括：

第一步：初始化参数，参数包括目标绩效阈值与数据绩效和阈值关系的大小关系判断方法、控制目标函数、剔除参数alpha_1、粘贴参数alpha_2和最小质量mass_min；将决策问题中根据因子变量{X1, X2, …, Xn}与多维目标{Y1, Y2, …, Yn}中某一个目标绩效指标的映射关系转换为抽象问题f(y|x)：根据因子变量集，预测满足某一个目标的绩效指标符合阈值要求的概率P(Y=1|X)；

第二步：求所有因子变量当前在alpha_1和1-alpha_1分位数的值，即统计当前数据点在各因子变量的累积分布概率等于alpha_1和1-alpha_1处的值，得到X（alpha_1）和X（1-alpha_1）；

第三步：选取各因子变量在分位数0至alpha_1、1-alpha_1至1之间的区域，生成备选的剔除空间区域；

第四步：判断是否有任一备选区域的剔除能使进行剔除后得到的空间区域的控制目标函数的结果比未进行剔除的原始空间区域的控制目标函数的结果大，是则进行下一步，否则跳转至第七步；

第五步：剔除能让目标函数结果增大最多的备选区域；

第六步：判断剩余区域形成新的超立方体的质量是否大于最小质量mass_min，是则进行下一步，否则跳转至第二步；

第七步：求所有因子变量当前向已剔除的区域往外扩展alpha_2分位数的值，即0-alpha_2和1+alpha_2分位数处的值，得到X（0-alpha_2）和X（1+alpha_2）；

第八步：选取各因子变量在分位数0-alpha_2至0、1至1+alpha_2之间的区域，生成备选的粘贴空间区域；

第九步：判断是否有任一备选区域的粘贴能使进行粘贴后得到的空间区域的控制目标函数的结果比未进行粘贴的原始空间区域的控制目标函数的结果大，是则粘贴让目标函数结果增大最多的备选区域并跳转至第七步，否则进行下一步；

第十步：输出案例数据库的目标规则集和因子变量关键变化范围；目标规则集包含控制目标函数值随各因子变量区域的剔除和粘贴的变化规律与特点，因子变量关键变化范围为最后剩余区域的变量上下限区间。

4.根据权利要求3所述的情景发现和脆弱性分析方法，其特征在于，所述第二步至第六步采用方法：

控制目标函数为：

；

其中，表示在m维空间的当前超立方体里的样本量，为案例数据库中的绩效结果，f(y)为当前超立方体里的正样本浓度，当前超立方体在各个维度上的上下限与样本群在该维度上的自变量上下限相同；

剔除参数alpha_1表示每次剔除的样本比例，低于分位数X（alpha_1）或高于分位数X（1-alpha_1）的样本将会被剔除；

依据剔除参数alpha_1，在当前超立方体内部生成一组备选的剔除空间区域：

其中：

；

；

和分别表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体里的α和1-α分位数的值；从所有的备选区域中选择一个，使得剔除该区域后，剩余空间区域的控制目标函数的结果比剔除前增加最多；从剔除一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选剔除区域为：

；

不断迭代执行剔除操作，直到备选区域集合中所有备选区域的剔除都不能增大PRIM控制目标函数结果，或满足以下条件：

；

指示函数的含义为：如果样本在当前超立方体里，则取值为1，反之为0，n是总样本量；是属于当前超立方体的样本量相对于总样本量的占比，即PRIM过程中超立方体的质量；是目标超立方体的最小样本占比，即预设的PRIM过程中超立方体的最小质量mass_min；当超立方体的质量达到了用户预先设定的最小质量时，PRIM的迭代剔除步骤结束。

5.根据权利要求4所述的情景发现和脆弱性分析方法，其特征在于，所述第七步至第十步采用方法：

根据粘贴参数alpha_2，在当前超立方体外部生成一组备选区域：

；

其中，

；

；

和表示因子变量{X1, X2, …, Xn}中的Xj在当前超立方体向外扩展α分位数的值，和表示Xj在中的最小值和最大值；从所有的备选区域中选择一个，使得粘贴该区域后，剩下区域的控制目标函数的结果比粘贴前增加最多；为粘贴一块区域后，可以得到一个新的超立方体：

；

其中，可选择的最优备选粘贴区域为：

；

不断迭代执行粘贴操作，直到备选区域集合中所有备选区域的粘贴都不能增大PRIM控制目标函数结果；迭代粘贴操作停止后，PRIM输出最终的超立方体作为结果，结果包含了案例数据库的目标规则集和因子变量的关键变化范围信息。

6.根据权利要求2-5任一所述的情景发现和脆弱性分析方法，其特征在于，所述帕累托最优集算法包括方法：

对多个目标采用以下公式进行优化：

；

；

；

；

其中o为多维目标F(x)={Y1, Y2, …, Yn}的目标个数，u为不等式约束的个数，v为等式约束的个数，和分别为中第i个变量x_i的上界和下界，w为变量x_i的总个数；

采用遗传算法，求帕累托最优集。

7.根据权利要求6所述的情景发现和脆弱性分析方法，其特征在于，所述遗传算法包括内容：

随机产生初始种群，采用非支配排序通过遗传算法的选择、交叉、变异三个基本操作得到第一代子代种群；

将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离迭代进行拥挤度计算：

；

为个体在目标空间中的曼哈顿距离；

对局部的支配关系排序为：

；

其中，和分别为局部空间中i和j两个个体的排序，和分别为i和j由曼哈顿距离所代表的个体拥挤度；

使用二元锦标赛选择方法或者参考点选择方法来维持种群的多样性；

根据非支配关系以及个体的拥挤度选取合适的个体组成新的父代种群；对于两个样本点，优先选择帕累托前沿排序更小的；如果两个样本点在同一个帕累托前沿上，则选择拥挤度更大的样本点；

判断是否满足预设的结束条件，是则结束，否则调整至将父代种群与子代种群合并，进行快速非支配排序，同时对每个非支配层中的个体按目标空间中的曼哈顿距离进行拥挤度计算步骤。

8.一种情景发现和脆弱性分析系统，应用于如权利要求1-7任一所述的情景发现和脆弱性分析方法，其特征在于，包括案例数据库、数据处理模块和数据输出模块；

所述案例数据库，存储对情景集的模拟产生的案例；

所述数据处理模块，判断案例数据库中的对情景集的模拟产生的案例的数量是否超过阈值；是则采用PRIM算法计算获得关键政策情景与脆弱性分析结果；否则采用帕累托最优集算法进行求解获得最优集后，判断最优集中案例数量是否超过阈值，是则对最优集采用PRIM方法计算获得关键政策情景与脆弱性分析结果，否则将最优集作为关键政策情景与脆弱性分析结果；

所述数据输出模块，将关键政策情景与脆弱性分析结果输出。

9.一种情景发现和脆弱性分析终端，包括存储器、处理器以及存储在所述存储器中并可在所述处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现如权利要求1至7任一项所述方法的步骤。

10.一种计算机可读存储介质，所述计算机可读存储介质存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现如权利要求1至7任一项所述方法的步骤。

Images