CN110390396B - 用于估计观测变量之间的因果关系的方法、装置和系统 - Google Patents

Abstract

公开了用于估计观测变量之间的因果关系的方法、装置和系统。根据本公开的方法，响应于接收到混合观测变量的观测数据，确定适用于连续观测变量和离散观测变量的混合因果关系目标式，其包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，且拟合不一致度基于观测变量的加权因子进行调整。然后利用混合观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对所述混合因果关系目标式进行最优化求解，以估计多个观测变量之间的因果关系。本公开的实施方式适用于混合观测变量的因果关系估计，并且因果网络结构对因观测变量估计误差而引起的敏感度较低，因此可以得到较为精确的因果关系。

Description

用于估计观测变量之间的因果关系的方法、装置和系统

技术领域

本公开涉及数据挖掘技术领域，更特别地涉及一种用于估计观测变量之间的因果关系的方法、装置和系统。

背景技术

在大数据时代，可以通过各种数据采集途径获得大量的数据。通过对这些数据进行数据分析与挖掘，可以从中得到很多有用的信息。然而，在众多应用领域中，人们往往只能看到系统的表象，却无法洞察系统背后复杂的作用机理及作用过程，而是仅能获得经验型的理解。

因果结构学习致力于基于系统的观测数据，自动还原系统背后复杂的作用机理，还原数据生成过程。目前，因果结构学习技术已被应用于制药、制造、市场分析等多个领域，以便洞察系统本质，继而指导决策，创造价值。

待被估计的观测变量可以具有不同数据类型的数据集，诸如包括连续变量的数据集或者离散变量的数据集。然而在现实的应用中，数据通常并非是由单一类型的变量组成，往往是连续变量和离散变量的混合。在目前的因果结构学习中，比较常用的模型包括结构方程模型、贝叶斯网络因果模型等，而这些常见的因果结构学习主要针对的是同一类型的数据集进行学习，对于混合数据集的情况其处理效果并不理想。

结构方程模型是当前用于因果结构学习的一种流行方法，例如LiNGAM模型，CAM模型，ANM模型等都属于结构方程模型。然而，这种类型的模型分析的是具有一致性的连续变量或者离散变量，并且目前的相关模型对于连续变量和混合变量的同时存在的情况处理效果很不理想。

贝叶斯网络因果模型是另一类被广泛运用的方法。传统的贝叶斯因果学习利用准确的推理算法支持因果结构的高精度求解。传统的贝叶斯网络因果模型在面对混合变量处理时通常采用连续变量离散化、基于条件分布概率的数据类型转换和基于条件独立性判定的方法。但是这些方法均存在一些问题。根据连续变量离散化的方法，将会把连续变量转换成离散变量，但是由于将连续变量离散化，这将损失连续变量的部分信息，因而存在精度问题。而根据基于条件概率分布的数据类型转换方法，将会把连续变量和离散变量的条件概率分布转换到同一种类型；但是由于要进行大量转换，因此时间复杂度高，无法支持观测变量维度较高的复杂因果结构学习任务。此外，对于基于条件独立性判定的方法，其将基于条件独立判定方法只能确定出网络结构的等价类，但是该方法并不能完全确定出网络中的具体因果关系。

此外，在Jing Xiang和Seyoung Kim等人发表的“A*Lasso for Learning aSparse Bayesian Network Structure for Continuous Variables”中，提出了一种利用A*Lasso算法来进行因果推理的方式。如图1所示，在该方案中，首先接收观测数据，然后基于使用利用观测变量进行拟合时的拟合不一致度和稀疏约束来执行稀疏因果结构建模，并利用A*Lasso搜索执行因果推理。此外，还通过设定候选因果序列数目的门槛值，来缩减最优因果序列搜索的搜索空间，来加速因果推理。

然而，目前已有的现有稀疏因果建模方法主要针对数据同为连续变量或者同为离散变量的数据集，而针对连续变量和离散变量共存的数据集，已有的因果模型或者推理算法要么时间复杂度高，无法支持观测变量维度较高的复杂因果结构学习，要么因连续变量离散化以及基于条件独立性判定而导致得到的结果并非最优结果，进而使得因果结构学习的精度下降。

为此，在本领域中需要一种针对混合变量情况来估计观测变量之间的因果关系的新技术方案。

发明内容

有鉴于此，本公开提出了一种用于估计观测变量之间的因果关系的方法、装置和系统，以至少部分上消除或者缓解现有技术中的问题。

根据本公开的第一方面，提供了一种用于估计观测变量之间的因果关系的方法。在该方法中，可以响应于接收到混合观测变量的观测数据，基于利用所述观测变量进行拟合时的拟合不一致度以及针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式。特别地，混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限，并且用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同。然后，可以利用所述混合观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系。

根据本公开的第二方面，提供了一种用于估计观测变量之间的因果关系的装置。该装置可以包括：混合目标式确定模块和混合因果推理模块。所述目标式确定模块可以被配置用于响应于接收到混合观测变量的观测数据，基于对利用观测变量进行拟合时的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式。特别地，所述混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限，并且用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同。所述混合因果推理模块可以被配置用于利用所述观测数据，通过分别用于所述离散类型观测数据的稀疏因果推理和用于所述连续类型观测数据的稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系。

根据本公开的第三方面，提供了一种用于估计观测变量之间的因果关系的系统，所述系统可以包括：处理器，和存储器，所述存储器存储了至少一个或多个计算机程序代码，所述计算机程序代码在被所述处理器执行时，使得所述处理器执行根据本公开的第一方面所述的方法的步骤。

根据本公开的第四方面，提供了一种计算机可读存储介质，其上包括有计算机程序代码，当由处理器执行中时，所述计算机程序代码使得所述计算机设备执行根据本公开的第一方面的方法的步骤。

根据本公开的第五方面，提供了一种计算程序产品，其上包括有计算机程序代码，当被加载到计算机设备中时，其可以使得该计算机设备执行根据本公开的第一方面的方法的步骤。

在本公开的实施方式中，针对混合观测变量的观测数据，建立既包括用于离散观测变量的因果关系目标式又包括用于连续观测变量的因果关系目标式的混合因果关系目标式，此外还通过观测变量的加权因子对拟合不一致度进行调整，并在推理过程中使用混合稀疏因果推理进行推理。因此，通过本公开的实施方式所提供的估计观测变量之间的因果关系的方案可以适用于混合观测数据情况下的因果关系估计，并且因果网络结构对因观测变量估计误差而引起的敏感度较低，因此可以得到更加精确的因果关系。而且，如将在下文中进一步的实施方式中所述，还可以采用基于变量组序列关系的搜索空间剪切技术来实现对搜索空间的剪切，以加速优化问题的求解。

附图说明

通过对结合附图所示出的实施方式进行详细说明，本公开的上述以及其他特征将更加明显，本公开的附图中相同的标号表示相同或相似的部件。在附图中：

图1示意性地示出了现有技术中的一种用于因果关系估计的技术方案；

图2示意性地示出了根据本公开的一个实施方式的用于估计观测变量之间的因果关系的方法的流程图；

图3示意性地示出了根据本公开的一个实施方式的用于稀疏因果关系推理的方法的流程图；

图4示意性地示出了根据本公开的另一实施方式的用于稀疏因果关系推理的方法的流程图；

图5示意性地示出了根据本公开的一个实施方式的用于生成变量组序列关系的方法的流程图；

图6A-6C示意性地示出了根据本公开的一个实施方式的生成变量组序列关系过程中的示例父关系图、强联通组件及其新有向无环图；

图7示意性地示出了根据本公开的一个实施方式的用于估计观测变量之间的因果关系的装置的流程图；以及

图8示意性地示出了根据本公开的一个实施方式的用于估计观测变量之间的因果关系的系统的示意图。

具体实施方式

在下文中，将参考附图详细描述本公开的各个示例性实施方式。应当注意，这些附图和描述涉及的仅仅是作为示例的优选实施方式。可以应该指出的是，根据随后的描述，很容易设想出此处公开的结构和方法的替换实施方式，并且可以在不脱离本公开要求保护的公开的原理的情况下使用这些替代实施方式。

应当理解，给出这些示例性实施方式仅仅是为了使本领域技术人员能够更好地理解进而实现本公开，而并非以任何方式限制本公开的范围。此外在附图中，出于说明的目的，将可选的步骤、模块、模块等以虚线框示出。

在此使用的术语“包括”、“包含”及类似术语应该被理解为是开放性的术语，即“包括/包含但不限于”。术语“基于”是“至少部分地基于”。术语“一个实施例”表示“至少一个实施例”；术语“另一实施例”表示“至少一个另外的实施例”。其他术语的相关定义将在下文描述中给出。

如前所述，现有的稀疏因果建模方法主要针对数据同为连续变量或者同为离散变量的数据集，而针对连续变量和离散变量共存的数据集，已有的因果模型或者推理算法要么时间复杂度高，无法支持观测变量维度较高的复杂因果结构学习，要么因连续变量离散化以及基于条件独立性判定而导致得到的结果并非最优结果，进而使得因果结构学习的精度下降。

为此，在本公开的实施方式中，提供了一种新的用于针对混合数据类型情况下的观测变量的因果关系估计的方法、装置和系统，以便至少部分上缓解上述问题。根据本公开的实施方式，提出针对混合变量的加权因果网络结构建模方法，其中综合观测数据的拟合不一致度、稀疏度约束和加权因子而建立混合因果网络结构，其中利用加权因子来消除不同数据类型的观测变量分析的敏感度和观测变量估计误差的敏感度。同时，基于新的混合因果模型，采用混合因果推理算法，来进行混合变量的因果结构学习。特别地，在一些实施方式中，可以从观测变量的可选父节点集合与当前子集的交集中选择最优父节点集合，而非从整个子集中进行选择，因而不仅能够保证搜索结果是最优结果，而且进一步加速推理过程。

在下文中，将参考图2至图8来描述根据本公开的用于估计观测变量的因果关系的方法、装置和系统。然而，需要说明的是，这些描述仅仅出于说明的目的，本公开并不局限于这些实施方式和附图中的细节。

图2示意性地示出了根据本公开的一个实施方式的用于估计观测变量之间的因果关系的流程图的示意图。如图2所示，首先在步骤201，响应于接收到混合观测变量的观测数据，基于利用所述观测变量进行拟合时的拟合不一致度以及针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式。其中所述混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式。所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限，并且用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同。

具体地说，可以设置观测数据库，在其中存储系统观测数据X，X∈R^N×D，X是一个N*D的矩阵，其N为观测样本的数量，D为观测变量的维数，或者观测变量的数目。观测数据库中的数据，可以是来自第三方的数据，或者是通过其他方式采集的数据。而且，可以预先对这些数据预处理，例如对原始数据的集成、规约、降噪等预处理。这些预处理操作本身在本领域是已知的，此处不再赘述。

在本公开的实施方式中，在构造因果关系目标式时，可以基于对利用观测变量进行拟合时的拟合不一致度以及针对所述观测变量的因果网络结构的稀疏度约束，来构造用于所述多个观测变量的混合因果关系目标式。但与传统的目标式不同的是，在本公开提出的方法中，将会建模一种适用于连续观测变量和离散观测变量的混合变量的因果模型，因果关系目标式包括用于离散类型观测数据的因果关系目标式和用于连续类型观测数据的因果关系目标式。而且将会利用观测变量的加权因子(也可简称为变量加权因子)对所述拟合不一致度进行调整，以便将不同类型的观测数据之间以及相同类型的观测数据之间的拟合不一致度的量级均调整到具有可比性的值范围内，从而降低因果网络结构对因观测变量估计误差而引起的敏感度。

在下文中，为了说明期间，将首先参考一个示例给出在本公开中提出的用于评价混合变量结构学习的拟合程度的混合信息条件准则(MIC)，该MIC是构建因果关系目标式的基础。例如，该MIC可以表示为：

在该式中的第一项是经过变量加权因子w_j调整的拟合不一致度，第二项是所设置的针对用来表征观测变量之间的因果关系的因果结构网络的稀疏度约束。具体地，G表示一个因果网络结构。x_i∈R^N是一个由N个观测值构成的向量，其描述了观测变量X_i；Pa(X_i)表示因果结构G中X_i的父节点，描述其他变量与X_i间的因果关系。w_i是一个标量参数，用于平衡变量间范围的差异，包括连续变量和离散变量的差异，也包括不同连续变量取值范围的差异；LL(X_i|Pa(X_i))基于数据的拟合不一致度，例如是负数对数似然值(negative log-likelihood)，也可称作损失函数)，Pen(X_i，Pa(X_i))用于设定对观测变量X_i因果估计时的稀疏约束，其例如可以是L0稀疏约束。

特别地，对于不同数据类型的观测变量，其拟合不一致度采用不同的方式来计算。对于连续变量，可以假设连续变量与其父节点是线性关系，且观测变量误差符合拉普拉斯分布。这样，则可以得到下式：

其中，β_i表示X_i和除X_i的其他变量X_-i之间的因果关系，其为待求解的因果关系,ε_i表示针对变量X_i的观测变量误差，Laplace(0,b_i)表示均值为0、尺度参数为bi的拉普拉斯分布。其中，如果β_ij＝0，则表示变量X_j不是X_i的父节点，因此(β_i，X_-i)即指示了X_i的父节点Pa(X_i)。

此处，假设变量误差分布是拉普拉斯分布，这是因为一方面超高斯分布对于奇异点而言更加鲁棒，另一方面在针对连续变量采用条件独立性判断的方法时，在误差并非高斯分布的情况下，可以得到可判定的因果结构。在拉普拉斯分布假设的情况下，针对混合变量中连续变量的MIC准则可以被表示为：

其中，n指示观测数据的数量，x_in指示变量X_i的第n个观测数据值，X_-i,n指示除变量X_i之外的其他变量X_-i的第n个观测数据值；bi指示针对变量X_i的拉普拉斯分布的尺度参数，β_i表示X_i和除X_i的其他变量X_-i之间的因果关系。

同样，对于离散变量，可以假设其误差服从逻辑(logistic)分布，进而可以得到下面的式子：

其中β_i表示X_i和除X_i的其他变量X_-i之间的因果关系，其为待求解的因果关系,ε_i表示针对变量X_i的变量误差，Logistic(0,1)表示均值为0、尺度参数为1的逻辑分布。在这种情况下，针对混合变量中连续变量的MIC准则可以被表示为：

其中n指示观测数据的数量，x_in指示变量X_i的第n个观测数据值，X_-i,n指示除变量X_i之外的其他变量X_-i的第n个观测数据值；β_i表示X_i和除X_i的其他变量X_-i之间的因果关系。

换句话说，在式子1中的混合信息条件判定准则MIC中，对于连续类型观测数据，将采用式(3)的MIC准则，而对于离散类型观测数据将采用式(5)的准则，而这两种准则的综合结果将作为混合变量结构学习的拟合程度的MIC准则。

此外，对于式(1)中的用来平衡变量间范围的差异的加权因子w_i，在下面给出了一种示例性的计算方法。

其中，Pa_c(X_i)表示对于X_i的最优可能父节点集合，X_i的父节点集合是最优可能父节点集合的一个子集，即因此，式1中的LL(X_i|Pa(X_i))/w_i是经过加权因此调整的相对损失，其相对于绝对损失而言，消减了数据不同变量范围引起的差异，同时削减模型对观测变量X_i估计误差的敏感度。

在本公开的实施方式中，针对不同类型的观测变量，上述式(6)中计算w_i所使用的具体式子是不同的。例如针对连续类型观测变量，可以拉普拉斯假设来计算。例如，w_i适于使用下式子来计算：

该式(7)与式(3)类似，但是不同在于使用了指示X_i的父节点的X’，其是X_-i的子集。类似地，对于离散类型观测变量，可以基于逻辑分布假设，例如w_i适于使用下式子来计算：

该式(7)与式(5)类似，但是不同在于使用了指示X_i的父节点的X’，其是X_-i的子集。

从式(7)和(8)还可以看出，所述变量加权因子的计算采用了与原有的目标式(无变量加权因子调整的目标式)类似的结构，用来计算在不考虑该变量xi的情况下，使用其他变量来拟合目标变量所需要的最小成本的下限，但是在计算变量加权因子的过程中，并不使用有向无环图约束，以便得到最小成本的下限。然而，需要说明的是，也可以采用其他类似的方式来确定变量加权因子，只要能够得到类似效果即可。

这样，基于上面的式子，可以针对每个观测变量Xi，首先通过上面的式子首先确定出其最优可能父节点集合，然后基于所确定的最优可能父节点集合而确定出相应地加权因子w_i。也就是说，其实在求解最优可能父节点集合情况下的负数对数似然值。因此，式子6也可以直接表示为：

进一步地，基于在本公开中提出的用于评价混合变量结构学习的拟合程度的混合信息条件准则(MIC)，可以给出下面的示例混合因果关系目标式，其中利用离散类型观测数据和连续类型观测数据的对数似然、变量加权因子和L0稀疏约束对稀疏因果结构进行鲁棒建模。

其中

x_i∈R^N是由N个观测值构成的向量，描述观测变量X_i；

x_-i∈R^N×(D-1)表示观测数据中除X_i变量外的所有数据，该矩阵中的每一列分别与变量X₁,…,X_i-1,X_i+1,…,X_D对应；

β_i∈R^D-1为待估计参数，描述变量X_i与其它变量间的因果关系强度；

||A||₀表示向量的L0范数；

约束G∈DAG用于约束因果网络结构是一个有向无环图；

λ_i为超参数，用于设定对观测变量X_j因果估计时的稀疏约束强度；

w_i为标量参数，为变量加权因子，用于削减模型对观测变量X_j

估计误差的敏感度；

z_i为数据类型指示，z_i＝1表示离散变量并且z_i＝0表示连续变量。

从上面的目标式(10)可以看出，该混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，针对离散观测变量，可以采用用于离散观测变量的因果关系目标式，并且对于连续观测变量，可以采用用于连续观测变量的因果关系目标式。此外，还利用变量分布用于连续观测变量和离散观测变量的加权因子w_i对所述拟合不一致度进行了调整，具体地将原始拟合不一致度与变量加权因子进行相除操作。然而需要说明的是，也可以采用其他类似的方式来进行调整，只要能够得到类似效果即可。

通过利用该变量加权因子w_i，可以将不同数据类型和相同数据类型的拟合不一致度转换为相互之间具有可比性的量级，从而消除目标式因为变量估计误差而引起的敏感度，使其更加鲁棒。

此外，从目标式(10)还可以看出，其中采用了L0稀疏约束。然而，需要说明的是，本发明并不局限于此，实际上，在目标式(1)中也可以采用L1稀疏约束。不过需要说明的是，L1稀疏约束在实现降维时，不仅会删除与目标不相关的因子，也有可能删除与目标相关的因子，从而影响稀疏降维的准确性，而L0稀疏约束则可以克服这一不足。

接着，在步骤202，利用所述混合观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系。

关于上述具有有向无环图约束的因果关系目标式，可以利用传统方法将其转化为一个最优因果序列递归求解问题。在根据本公开的一个实施方式中，可以将转换后的递归问题描述如下：

其中，V描述所有观测变量的集合，U描述尚未排序观测变量的集合，描述隶属于变量集合U且在因果序列中排序在第一个位置的变量。S(β_j)表示向量β_j的支撑集合。

针对式(11)的求解存在各种方法，比如传统的A*搜索算法或者任何其他适当的算法。但是对于连续数据和离散数据分布采用适合他们各自的特征选择标准。

此外，在求解所述的因果序列递归问题时，也可以使用扩展的A*搜索算法，其为针对传统的A*搜索算法的一种扩展。在被公开的一个实施方式中，提出了基于一种扩展A*搜索算法来进行求解。该扩展A*搜索算法在此处可称为混合A*FoBa算法，用于在因果序列图空间中搜索最优的因果序列生成路径。但是对于zi＝1的情况，离散变量的情况，将采用离散变量相关的FoBa稀疏回归算法，比如在特征选择时基于F1-score；而对于zi＝0的情况，即连续变量的情况，将采用连续变量相关的FoBa稀疏回归算法。此时，在特征选择时比如可以基于平均绝对值误差MAE。

在下文中，将根据该混合A*FoBa算法，f函数可以表示如下：

f(Q_s)＝g(Q_s)+h(Q_s) (式13)

其中，g(Q_s)表示从起始节点开始经过该节点Q_s至当前节点所花费的成本，而h(Q_s)表示从当前节点至目标节点的预计成本，其中g(Q_s)和h(Q_s)可以分别按照下式计算，

需要说明的是，为保证扩展A*搜索算法能找到最优解，h()函数应该满足两个特性：一致性(consistent)及可接受性(admissible)。其中，一致性是假设对任意两个节点Q_s和Q_s’，Q_s’＝[Q_s,x_j],如果h(Q_s)<＝h(Q_s’)+MScore(x_j|Q_s)成立，则称h()函数满足一致性。可接受性是指，对从当前节点Q_s到目标节点所需成本的估计h(Q_s)，应该是对真实成本的欠估计。

从式(6)可以看出在混合A*FoBa搜索算法中，h函数和g函数与传统的A*搜索算法不同，他们基于经过变量加权因子调整的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，来确定从起始节点至当前节点的所花费成本(g函数)和从当前节点至目标节点的预测成本(h函数)。

在下文中，将参考图3来描述根据公开中的中的混合A*FoBa搜索方法的示例流程。

参考图3，首先在步骤301执行初始化操作。具体地，针对因果序列集合CLOSE、已排序变量集合Q_s、候选因果序列和因果序列集合OPEN进行初始化。

因果序列集合CLOSE用于存放已经根据其生成新的候选排序因果序列的已排序序列集合，其可以在初始化时被初始化为空集合。已排序变量集合Q_s是已被排序的变量的集合，在初始化时可以被初始化为空集合。候选因果序列为已经排序的变量的候选顺序。因果序列集合OPEN中可以采用下面的形式：

其中Q_s是已排序变量集合,f(Q_s)是上面的式(13)，是候选因果序列，其可以被初始化为：

接着，在步骤302，从OPEN集合中抽取f函数值最小的已排序变量集合Qs，并在步骤303确定Qs的h函数是否为0，即是否已经找到一条到达目标结点的路径。如果为零，则将与f函数最小值对应的候选因果序列作为最优因果序列返回，否则方法进入步骤304。在初始执行时，在OPEN集合中并不存在任何已排序变量集合，此时方法也将进行至步骤304，以生成新的已排序变量集合。

在步骤304，针对每个观测变量v∈V\Qs，生成新的已排序因果关系序列集合Q_s’和候选因果序列并学习对应的f(Q_s’),

对于初始执行，Q_s为空的情况下，可以将针对每个观测变量v生成新的已排序因果集合、相应的候选因果序列并学习对应的f(Q_s’)。

然后在步骤305，将生成的候选因果序列添加到OPEN集合中，以更新OPEN集合，即:

接下来，可以在步骤306，将之前抽取的f函数值最小的Q_s从OPEN集合中转移到CLOSE集合中，即将Q_s及其相关信息从OPEN集合中删除，并添加到CLOSE集合中。然后，该方法跳转到302继续执行。

步骤301至305的操作可以反复执行，直至找到h函数为0的已排序变量集合，并将其对应的候选因果序列作为最优因果序列返回。

出于说明的目的，将以三个观测变量为例进行简要说明。然而需要说明的是实际应用中，观测变量的维度会高得多。例如对于三个观测变量v₁,v₂和v₃，在初始执行OPEN集合为空，此时将针对每个观测变量v₁,v₂和v₃，生成新的已排序观测变量集合{v₁},{v₂}和{v₃}，获得相应的候选因果序列并学习对应的f(Q_s’)。生成的候选因果序列将会被放在OPEN集合中。在第二次执行的时候，从这三个变量集合中找出f函数值最小的已排序变量集合，例如{v₂}，然后判断其h函数值，此时其h函数值并不为零，因此接着生成例如{v₁，v₂}和{v₂，v₃}两个已排序集合及其对应的候选因果排序v₂，v₁和v₂，v₃和对应的f(Q_s’)。然后{v₂}及其相关信息从OPEN集合中删除，并放入CLOSE集合中。在第三次执行时，从OPEN集合中继续选择f函数最小的序列，比如是{v₂，v₃}，进一步生成已排序集合{v₁，v₂，v₃}及其对应的候选因果排序v₂，v₃，v₁和对应的f(Q_s’)，并将{v₂，v₃}及其相关信息从OPEN集合中转入CLOSE集合中。在第四次执行时，再次选出选择f函数最小的序列，此时是{v₁，v₂，v₃}，且其h函数值为零，因此可以将对应的候选因果序列v₂，v₃，v₁作为最优因果序列返回，从而完成搜索过程。

此外，在进行迭代操作的搜索过程中还可以采用搜索空间剪切来加速稀疏因果关系推理的速度。与现有技术中的基于门槛值的搜索空间剪切不同，在本公开中的一个实施方式中，采用了基于变量组序列关系的搜索空间剪切。换句话说，在本公开中提议引入变量组序列关系并以该变量组序列关系来对序列图搜索空间进行剪切。图4示出了在采用了搜索空间剪切情况下的混合A*FoBa搜索方法的示例流程图。

参考图4，图4中的步骤401至407与图3中的步骤301至307是相同的，此处不再赘述；其不同之处在于在步骤404和405之间，新增了一个基于变量组序列关系的搜索空间剪切步骤408。当在步骤404得到新的已排序因果关系序列集合Q_s’和对应的候选因果序列后将在步骤408进一步利用变量组序列关系对该候候选因果序列进行核验，以确定其是否与变量组序列关系相冲突。如果该候选因果序列与变量组序列关系不存在冲突，则步骤进行至步骤405，将新生成的候选因果序列加入OPEN集合。相反，如果该候选因果序列与变量组序列关系存在冲突，则将跳过步骤405直接执行步骤406。通过此方式，可以对搜索空间进行剪切，以便加速最优求解过程。

该变量组序列关系可以通过各种方式获取。例如，可以由相关领域专家给出变量组序列关系，例如对于G₁＝{v₁,v₂},G₂＝{v₃,v₆},G₃＝{v₃,v₅,v₇,v₈}，领域专家可以给出变量组序列关系另外，也可以通过对观测数据的自动学习来得到变量组关系。在下面，将参考图5来描述一种自动学习变量组序列关系的一种示例性方法。

参考图5，首先在步骤501，针对每个观测变量v_i(i∈{1，...，D})，获取其潜在最优父结点集合POP_i。可以使用各种方法来获取各个观测变量的潜在最优父节点集合。例如可以使用传统的A*搜索稀疏回归算法、A*lasso搜索算法等。此外，也可以采用本公开中提供的混合A*FoBa稀疏回归算法获取该集合。获取最优父节点集合的具体方式在本领域是已知的，此处不再赘述。

接着，在步骤502，基于在步骤501获取的最优父结点集合{POPS₁，...，POPS_D},来生成父关系图。具体地，可以将每个观测变量转换成图中的一个节点，针对每个观测变量V_p∈POPS_i,可以生成一条由V_p指向V_i的边，以描述变量V_p是变量V_i的一个潜在父节点。图6A示出了根据本公开的实施方式的父关系图的一个示例，其中涉及6个变量。

然后，在步骤503，从步骤502生成的父关系图中抽取强联通组件SSC，得到{SSC₁，...，SSC_M}，其中SSC₁∪...∪SSC_M＝{v₁，...，v_D}，即所有强联通组件的并集是所述变量的全集。强联通组件是自身能够形成一个环形回路且与其它节点或强联通组件不能形成环形回路的节点的集合。图6B示意性地示出了从图6A所示父关系图抽取的强联通组件SCC₁和SCC₂。

接下来，在步骤504，基于生成的强联通组件SSCs，将步骤502父关系图转换为一个新的有向无环图G_new。一个示例方法是将每个强联通组件SSC_i转换成一个新结点N_i，如果在父关系图中存在SSC_i中的任意节点指向SSC_j中的节点的边，则在G_new中添加N_i指向N_j的边。图6C示出了根据图6B中的强联通组件所生成的新有向无环图。图6B中的两个强联通组件SCC₁和SCC₂分别被转换为两个新节点N₁和N₂。由于在图6A的父关系图中，在SCC₁中存在指向SCC₂中的边，所以在两个节点N₁和N₂之间增加一个N₁指向N₂的边。以此方式，可以得到新的有向无环图G_new。

最后，在步骤505，抽取所述有向无环图的强联通组件，将其序列关系作为变量组序列关系。可以直接有向无环图中的序列关系从而获取变量组序列关系例如对于图6C所示的新有向无环图，可以得到

以此方式，就可以得到用于进行搜索空间剪切的变量组序列关系。然而，需要说明的是，本发明并不局限于此，任何其他能够得到变量组序列关系的方式也是可行的。

在本公开的实施方式中，针对混合观测变量的观测数据，建立既包括用于离散观测变量的因果关系目标式又包括用于连续观测变量的因果关系目标式的混合因果关系目标式，此外还通过观测变量的加权因子对拟合不一致度进行调整，并在推理过程中使用混合稀疏因果推理进行推理。因此，通过本公开的实施方式所提供的估计观测变量之间的因果关系的方案可以适用于混合观测数据情况下的因果关系估计，并且因果网络结构对因观测变量估计误差而引起的敏感度较低，因此可以得到更加精确的因果关系。而且，如将在下文中描述的进一步的实施方式中所述，还可以采用基于变量组序列关系的搜索空间剪切技术来实现对搜索空间的剪切，加速优化问题的求解。。

图7还示出了根据本公开的用于估计观测变量之间的因果关系的示例装置700的方框图。如图7所示，装置700包括混合目标式确定模块702和混合因果推理模块703。

所述目标式确定模块702可以被配置用于响应于接收到混合观测变量的观测数据，基于对利用观测变量进行拟合时的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式。所述混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式。所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限。并且，用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同。可以设置观测数据库701，在其中存储系统观测数据X，X∈R^N×D，X是一个N*D的矩阵，其N为观测样本的数量，D为观测变量的维数，或者观测变量的数目。观测数据库中的数据，可以是来自第三方的数据，或者是通过其他方式采集的数据。而且，可以是预先经过处理的数据。

所述因果推理模块703被配置为被配置为利用所述混合观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系。

在根据本公开的实施方式中，用于连续观测变量的因果关系目标式可以基于观测变量误差符合拉普拉斯分布的假设而建立，并且用于离散观测变量的因果关系目标式可以基于观测变量误差符合逻辑分布的假设而建立。

再根据本公开的实施方式中，用于连续观测变量的加权因子基于观测变量误差符合拉普拉斯分布的假设而计算，并且用于离散观测变量的加权因子基于观测变量误差符合逻辑分布的假设而建立。

在根据本公开的实施方式中，所述混合因果推理模块被配置为将针对所述多个观测变量的变量模型的混合稀疏因果推理转换为最优因果序列递归求解的问题。所述因果推理模块703可以被配置为基于经过调整的所述拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，来确定从起始节点至当前节点的所花费成本(g函数值)和从当前节点至目标节点的预测成本(h函数值)。关于最优因果序列递归求解的问题，可以参考上文中参考图3的描述。

在根据本公开的一个实施方式中，所述稀疏度约束可以是L0稀疏约束。

在根据本公开的一个实施方式中，所述装置700还可以进一步包括搜索空间剪切模块704，其被配置为在所述最优因果序列递归求解的过程中针对生成新的候选因果序列，确定其是否与预定的变量组序列关系相冲突；以及如果确定存在冲突，则丢弃所述生成的新候选因果序列。关于搜索空间剪切，可以参考上文中参考图4的描述。

在根据本公开的一个实施方式中，所述预定的变量组序列关系是由专家给定的。在本公开的另一实施方式中，所述预定的变量组序列关系是基于观测数据自动确定的。

在根据本公开的另一实施方式中，所述装置700还可以进一步包括序列关系确定模块705。该序列关系确定模块705可以被配置为通过下述操作自动确定所述预定的变量组序列关系：针对所述观测变量中的每个观测变量，获取其潜在最优父节点集合；基于获取的所述最优父节点集合，生成父关系图；从所述父关系图中抽取强联通组件；通过将所述强联通组件中的每个组件转换成新节点，并当两个强联通组件在父关系图中具有联通性时，在对应的新节点之间增加相应的边，来将所述父关系图转换成新的有向无环图；以及提取所述有向无环图中的强联通组件之间的序列关系，以作为所述预定的变量组序列关系。关于自动变量组序列关系的确定，可以参考上文参照图5的描述。

图8示意性地示出了根据本公开的一个实施方式的用于估计观测变量之间的因果关系的系统的示意图。下面，将参考图8来描述可以实现本发明的因果关系估计的系统。

图8中所示的计算机系统包括CPU(中央处理单元)801、RAM(随机存取存储器)802、ROM(只读存储器)803、系统总线804、硬盘控制器805、键盘控制器806、串行接口控制器807、并行接口控制器808、显示器控制器809、硬盘810、键盘811、串行外部设备812、并行外部设备813和显示器814。在这些部件中，与系统总线804相连的有CPU 801、RAM 802、ROM 803、硬盘控制器805、键盘控制器806、串行接口控制器807、并行接口控制器808和显示器控制器809。硬盘810与硬盘控制器805相连，键盘811与键盘控制器806相连，串行外部设备812与串行接口控制器807相连，并行外部设备813与并行接口控制器808相连，以及显示器814与显示器控制器809相连。

在存储器中可以存储一个或多个代码，所述代码在被所述计算机执行时，指示所述CPU执行在本公开的实施方式中提出的方法的步骤，例如前面参考图2至图6所描述的那些是实施方式。

需要说明的是，图8所述的结构方框图仅仅为了示例的目的而示出的，并非是对本发明的限制。在一些情况下，可以根据需要添加或者减少其中的一些设备。

可以理解的是，在本公开提出的方案可以在制药、制造、市场分析、交通预测、天气预测、空气质量预测等各种应用中使用，以获得有益效果。

此外，本发明的实施方式可以以软件、硬件或者软件和硬件的结合来实现。硬件部分可以利用专用逻辑来实现；软件部分可以存储在存储器中，由适当的指令执行系统，例如微处理器或者专用设计硬件来执行。

本领域的普通技术人员可以理解上述的方法和装置可以使用计算机可执行指令和/或包含在处理器控制代码中来实现，例如在诸如磁盘、CD或DVD－ROM的载体介质、诸如只读存储器(固件)的可编程的存储器或者诸如光学或电子信号载体的数据载体上提供这样的代码。

本实施方式的装置及其组件可以由诸如超大规模集成电路或门阵列、诸如逻辑芯片、晶体管等的半导体、或者诸如现场可编程门阵列、可编程逻辑设备等的可编程硬件设备的硬件电路实现，也可以用由各种类型的处理器执行的软件实现，也可以由上述硬件电路和软件的结合例如固件来实现。

虽然已经参考目前考虑到的实施方式描述了本发明，但应该理解本发明不限于所公开的实施方式。相反，本发明旨在涵盖所附权利要求书的精神和范围内所包括的各种修改和等同布置。所附权利要求书的范围符合最广泛解释，以包含所有这样的修改及等同结构和功能。

Claims (18)

1.一种由计算机执行的用于交通预测的方法，包括：

将混合观测变量的交通观测数据采集到观测数据库中，其中所述交通观测数据与交通有关并且用于所述交通预测；

从所述观测数据库中接收所述交通观测数据，所述交通观测数据包括多个观测样本，所述混合观测变量包括与交通有关的连续观测变量和与交通有关的离散观测变量，所述连续观测变量和/或所述离散观测变量指示对应观测样本的一组观测值；

基于利用所述观测变量进行拟合时的拟合不一致度以及针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式，其中所述混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，其中所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限，并且用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同；

利用所述交通观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系；以及

基于所述交通观测数据以及所述因果关系进行交通预测。

2.根据权利要求1所述的方法，其中用于连续观测变量的因果关系目标式基于观测变量误差符合拉普拉斯分布的假设而建立，并且用于离散观测变量的因果关系目标式基于观测变量误差符合逻辑分布的假设而建立。

3.根据权利要求1所述的方法，其中用于连续观测变量的加权因子基于观测变量误差符合拉普拉斯分布的假设而计算，并且用于离散观测变量的加权因子基于观测变量误差符合逻辑分布的假设而建立。

4.根据权利要求1至3任一项所述的方法，其中将针对所述多个观测变量的混合因果关系目标式的混合稀疏因果推理转换为最优因果序列递归求解的问题，其中基于经过所述观测变量的加权因子调整的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，来确定在已排序因果序列中从起始节点至当前节点所花费的成本和从当前节点至目标节点的预测成本。

5.根据权利要求1至3任一项所述的方法，其中所述稀疏度约束是L0稀疏约束。

6.根据权利要求4所述的方法，进一步包括：在所述最优因果序列递归求解的过程中：

针对生成新的候选因果序列，确定其是否与预定的变量组序列关系相冲突；以及

如果确定存在冲突，则丢弃生成的所述新的候选因果序列。

7.根据权利要求6所述的方法，其中所述预定的变量组序列关系基于以下任何一个或多个方式来确定：

由领域专家给定；以及

基于所述交通观测数据自动确定。

8.根据权利要求7所述的方法，其中所述预定的变量组序列关系通过下述操作而基于所述交通观测数据自动确定：

针对所述观测变量中的每个观测变量，获取其潜在的最优父结点集合；

基于获取的所述最优父结点集合，生成父关系图；

从所述父关系图中抽取强联通组件；

通过将所述强联通组件中的每个组件转换成新节点，并当两个强联通组件在父关系图中具有联通性时，在对应的新节点之间增加相应的边，来将所述父关系图转换成新的有向无环图；以及

提取所述有向无环图中的强联通组件之间的序列关系，以作为所述预定的变量组序列关系。

9.一种用于交通预测的装置，包括：

混合目标式确定模块，被配置为从其中采集了混合观测变量的交通观测数据的观测数据库中接收所述交通观测数据，其中所述交通观测数据与交通有关并且用于所述交通预测，所述交通观测数据包括多个观测样本，所述混合观测变量包括与交通有关的连续观测变量和与交通有关的离散观测变量，所述连续观测变量和/或所述离散观测变量指示对应观测样本的一组观测值；基于利用所述观测变量进行拟合时的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，确定适用于连续观测变量和离散观测变量的混合因果关系目标式，其中所述混合因果关系目标式包括用于连续观测变量的因果关系目标式和用于离散观测变量的因果关系目标式，其中所述拟合不一致度基于所述观测变量的加权因子进行调整，所述观测变量的加权因子表示使用除该观测变量之外的其他观测变量来拟合目标变量所需的最小成本的下限，并且用于所述离散观测变量的加权因子与用于所述连续观测变量的加权因子不同；以及

混合因果推理模块，被配置为利用所述混合观测数据，通过适用于连续观测变量和离散观测变量两者的混合稀疏因果推理，在有向无环图约束下，对用于所述多个观测变量的混合因果关系目标式进行最优化求解，以估计所述多个观测变量之间的因果关系，以及基于所述交通观测数据以及因果关系进行交通预测。

10.根据权利要求9所述的装置，其中用于连续观测变量的因果关系目标式基于观测变量误差符合拉普拉斯分布的假设而建立，并且用于离散观测变量的因果关系目标式基于观测变量误差符合逻辑分布的假设而建立。

11.根据权利要求9所述的装置，其中用于连续观测变量的加权因子基于观测变量误差符合拉普拉斯分布的假设而计算，并且用于离散观测变量的加权因子基于观测变量误差符合逻辑分布的假设而建立。

12.根据权利要求9至11其中任一项所述的装置，其中所述混合因果推理模块被配置为将针对所述多个观测变量的混合因果关系目标式的混合稀疏因果推理转换为最优因果序列递归求解的问题，其中所述因果推理模块被配置为基于经过所述观测变量的加权因子调整的拟合不一致度和针对所述观测变量的因果网络结构的稀疏度约束，来确定在已排序因果序列中从起始节点至当前节点所花费的成本和从当前节点至目标节点的预测成本。

13.根据权利要求9至11其中任一项所述的装置，其中所述稀疏度约束是L0稀疏约束。

14.根据权利要求12所述的装置，进一步包括搜索空间剪切模块，其被配置为在所述最优因果序列递归求解的过程中：

针对生成新的候选因果序列，确定其是否与预定的变量组序列关系相冲突；以及

如果确定存在冲突，则丢弃生成的所述新的候选因果序列。

15.根据权利要求14所述的装置，其中所述预定的变量组序列关系基于以下任何一个或多个方式来确定：

由领域专家给定；以及

基于所述交通观测数据自动确定。

16.根据权利要求15所述的装置，还进一步包括：序列关系确定模块，其被配置为通过下述操作自动确定所述预定的变量组序列关系：

针对所述观测变量中的每个观测变量，获取其潜在的最优父节点集合；

基于获取的所述最优父节点集合，生成父关系图；

从所述父关系图中抽取强联通组件；

通过将所述强联通组件中的每个组件转换成新节点，并当两个强联通组件在父关系图中具有联通性时，在对应的新节点之间增加相应的边，来将所述父关系图转换成新的有向无环图；以及

提取所述有向无环图中的强联通组件之间的序列关系，以作为所述预定的变量组序列关系。

17.一种用于交通预测的系统，包括：

处理器；以及

存储器，其中存储有计算机程序代码，所述计算机程序代码在被所述处理器执行时，使得所述处理器执行根据权利要求1至8其中任一项所述的方法。

18.一种计算机可读存储介质，其上存储有计算机程序代码，所述计算机程序代码在被处理器执行时，执行根据权利要求1至8其中任一项所述的方法。