具体实施方式

本发明的核心是提供一种级联医疗观测数据的因果推断方法及系统，能够很好地识别具有级联结构的医疗观测数据的因果方向，显著提高因果方向识别的精确度。

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

本发明实施例一方面提供一种级联医疗观测数据的因果推断方法，请参阅图1，该方法包括：

步骤S110、获取级联医疗观测数据，并从级联医疗观测数据中抽取第一变量和第二变量。

本发明实施例中，在获取到的级联医疗观测数据中，包括观测到的作为初始原因和作为最终结果的数据，另外在初始原因和最终结果之间还存在未被观测到的中间变量和加性噪声，我们只观察作为初始原因和最终结果的数据，提取出来作为第一变量和第二变量。其中，第一变量和第二变量具有因果关系，且正确的因果方向不确定，既可以是第一变量-第二变量，即第一变量作为原因、第二变量作为结果；也可以是第二变量-第一变量，即第二变量作为原因、第一变量作为结果。

步骤S120、以因果关系中的原因、级联结构中每个深度对应的中间变量以及因果关系中的结果为参数，建立改进后的级联非线性加性噪声模型。

需要说明的是，ANM(英文全称：Additive Noise Model，中文全称：加性噪声模型)模型是一种非线性情况下二元变量之间常用的因果发现算法，其模型可形式化为y＝f(x)+ε，ε⊥x。而CANM(英文全称：Cascade Additive Noise Model，中文全称：级联非线性加性噪声模型)是一种为了研究间接和非线性变量之间的因果关系而提出的模型，主要是将ANM模型应用在含有中间变量的数据上。

本发明实施例中，改进后的级联非线性加性噪声模型可以看做是由多个ANM模型的组合，每个直接的因果影响都遵循ANM模型，且因果关系中的原因和结果之间存在未被观测到的中间变量和潜在噪声，该模型可以更好地匹配具有级联结构的医疗观测数据。

步骤S130、针对改进后的级联非线性加性噪声模型，构造边缘对数似然函数。

步骤S140、对边缘对数似然函数进行对抗变分解，并利用近似后验分布方法优化得到边缘对数似然函数对应的变分下界。

本发明实施例中，对边缘对数似然函数进行对抗变分解，并利用简单分布来近视难以求解的真实后验分布，通过对抗的方法优化得到边缘对数似然函数对应的变分下界，将最大化边缘对数似然转化为最大化变分下界。

步骤S150、以第一变量作为因果关系中的原因，以第二变量作为因果关系中的结果，利用预置的对抗训练模型求解最大化的变分下界，得到第一变分下界值；以第二变量作为因果关系中的原因，以第一变量为作为因果关系中的结果，利用对抗训练模型求解最大化的变分下界，得到第二变分下界值。

本发明实施例中，首先假定级联医疗观测数据中正确的因果方向是第一变量-第二变量，则可以将第一变量作为因果关系中的原因，第二变量作为结果，通过预置的对抗训练模型求解最大化的边缘对数似然函数对应的变分下界，得到第一变分下界值；然后假定级联医疗观测数据中正确的因果方向是第二变量-第一变量，则相应地将第二变量作为因果关系中的原因，第一变量作为结果，通过同样的方法可以得到第二变分下界值。其中，在通过预置的对抗训练模型求解最大化的边缘对数似然函数对应的变分下界时，我们利用对抗策略而非近似公式来绕开KL散度，能够允许加性噪声有更广泛的分布，从而提升模型的推断能力。

步骤S160、将第一变分下界值与第二变分下界值进行比较，得到比较结果，并根据比较结果，确定级联医疗观测数据的因果方向。

本发明实施例中，求解得到第一变分下界值与第二变分下界值后，将两个值进行比较，若第一变分下界值大于第二变分下界值，则确定级联医疗观测数据中正确的因果方向是第一变量-第二变量；反之，则确定级联医疗观测数据中正确的因果方向是第二变量-第一变量。

以上可知，本发明实施例提供的级联医疗观测数据的因果推断方法，以因果关系中的原因、级联结构中每个深度对应的中间变量以及因果关系中的结果为参数，建立改进后的级联非线性加性噪声模型，可以更好地匹配具有级联结构的医疗观测数据，提高识别级联医疗数据因果方向的精确度，同时，通过预置的对抗训练模型求解最大化的边缘对数似然函数对应的变分下界，利用对抗策略而非近似公式来绕开KL散度，能够允许加性噪声有更广泛的分布，从而提升模型的推断能力，与现有技术相比，能够很好地识别具有级联结构的医疗观测数据的因果方向，显著提高因果方向识别的精确度。

请参阅图2，作为本发明一种优选实施例，改进后的级联非线性加性噪声模型的表达式为：

Z₁＝f₁(X；θ)+ε₁

Z_T＝f_T(Z_pa(r)；θ)+ε_r

Y＝f_T+1(Z_pa(y)；θ)+ε_y

式中，T代表级联结构的深度，X代表因果关系中的原因，Z_T代表级联结构中每个深度对应的中间变量，Y代表因果关系中的结果，f＝{f₁，f₂，.....，f_T}代表一组非线性函数，θ代表因果关系中的参数，ε_T代表级联结构中每个深度对应的加性噪声，Z_pa(T)代表在级联结构中Z_T的前一个深度对应的中间变量，Z_pa(y)代表在级联结构中最后一个深度对应的中间变量，ε_y代表从Z_pa(y)到Y的加性噪声。

本发明实施例中，假设因果机制中没有混杂因素，且数据生成过程中遵循非线性加性噪声假设，则因果关系中的原因X、级联结构中每个深度对应的加性噪声ε_T和从Z_pa(y)到Y的加性噪声ε_y相互独立。

进一步地，上述实施例中，边缘对数似然函数的表达式为：

式中，p_θ()代表似然函数，xⁱ代表X中的第i个数据点，yⁱ代表Y中的第i个数据点，z代表中间变量，其中，i＝1，2，3...m，m代表数据点的个数。

本发明实施例中，假设数据为原因X包括m个数据点xⁱ，原因Y包括m个数据点yⁱ，则可以得到数据D的边缘对数似然函数的表达式。

更进一步地，上述实施例中，步骤S140包括：

利用马尔科夫条件对边缘对数似然函数进行分解，得到分解后的边缘对数似然函数的表达式：

分别对上述表达式中的p_θ(yⁱ|z_pa(y))和p_θ(z_t|z_pa(t))进行分解，并将函数f_T+1(Z_pa(y))改写为f(x，ε)，得到改写后的边缘对数似然函数的表达式：

式中，代表结果变量的加性噪声；ε代表中间变量的加性噪声。

本发明实施例中，首先利用马尔科夫条件对边缘对数似然函数进行分解，然后通过原因与加性噪声之间的独立性分解p_θ(yⁱ|z_pa(y))和p_θ(z_t|z_pa(t))，同时，由于最后一个未被观测到的中间变量Z_pa(y)包含了所有加性噪声ε_T和原因X对结果Y的所有影响，所以可以将函数f_T+1(Z_pa(y))改写为f(x，ε)，之后可以得到改写后的边缘对数似然函数的表达式。

引入参数利用简单分布来近似后验分布p_θ(ε|xⁱ，yⁱ)，对所述边缘对数似然函数进行进一步分解，得到进一步分解后的所述边缘对数似然函数的表达式：

将上述表达式中的第一项定义为变分下界，则当时，上述表达式中的KL散度为0，边缘对数似然函数等于边缘对数似然函数对应的变分下界，对边缘对数似然函数对应的变分下界进行分解，得到分解后的边缘对数似然函数对应的变分下界的表达式：

本发明实施例中，利用一个关于参数形式简单的分布来近似关于参数θ而难以求解的真实后验分布p_θ(ε|xⁱ，yⁱ)，通过对抗的方法利用和p_θ(ε|xⁱ，y^I)联合优化数据D的边缘对数似然对应的变分下界(ELBO)，且边缘对数似然为每个数据点(xⁱ，yⁱ)的边缘对数似然之和。由于我们希望能够更好地近似出p_θ(ε|xⁱ，yⁱ)，所以要求最小，但是对于数据D的边缘对数似然之和已经固定，所以我们需要最大化这一项称为变分下界，且当时，边缘对数似然函数等于边缘对数似然函数对应的变分下界，所以最大化边缘对数似然等价于最大化边缘对数似然对应的变分下界。

将上述表达式中的最后一项改写为构造判别网络模型T(X，Y；ε)，隐式的将表述为判别网络模型T(X，Y；ε)的最优值，利用判别网络的对抗策略来绕过KL散度，进一步得到所述边缘对数似然函数对应的变分下界的表达式：

式中，T^*(X，Y；ε)代表判别网络模型T(X，Y；ε)的最优值。

本发明实施例中，由于计算KL散度除了少数分布外，其积分项不具有封闭形式的解析解，所以我们将上述表达式中的最后一项改写为并通过隐式的将表述为我们构造的一个判别网络模型T(X，Y；ε)的最优值，利用判别网络的对抗策略来绕过KL散度，对抗策略能够允许更广泛的分布作为潜在噪声的先验，并使数据D到加性噪声的映射更灵活。所以，我们可以进一步得到边缘对数似然函数对应的变分下界的表达式。

请参阅图3，作为本发明一种优选实施例，对抗训练模型采用带有判别网络的变分自动编码器，包括编码器模块、解码器模块和判别器模块。

本发明实施例中，对抗训练模型由编码器模块、解码器模块和判别器模块这三部分组成，通过使用编码器模块、解码器模块和判别器模块的交替处理，来优化边缘对数似然函数对应的变分下界使其收敛，求解得到第一变分下界值和第二变分下界值，以确定带有中间变量的因果数据之间的因果方向。

进一步地，上述实施例中，编码器模块将简单分布表述为编码网络，编码网络采用三个带有ReLU非线性函数的全连接层和一个不做非线性处理的输出层作为网络结构，通过编码网络将因果关系中的原因和预置的随机变量映射为加性噪声。

本发明实施例中，编码器模块将简单分布表述为编码网络，即encoder，由于编码网络encoder是医疗观测级联数据到加性噪声ε的映射模型，所以我们使用三个带有ReLU非线性函数的全连接层和一个不做非线性处理的输出层作为编码网络encoder的网络结构。通过编码网络encoder将医疗观测级联数据和随机变量u(u～N(0，l))一起编码为加性噪声ε，而不使用重参数化，这样能够使医疗观测级联数据到加性噪声的映射更加灵活，能让模型学习更复杂的概率分布。

更进一步地，上述实施例中，解码器模块将后验分布p_θ(y′ⁱ|xⁱ，ε)表述为解码网络，解码网络采用与编码网络相同的网络结构，通过解码网络将因果关系中的原因和编码网络输出的加性噪声重构为结果，并计算重构后的结果与因果关系中的结果之间的重构误差，采用蒙特卡罗方法估计变分下界的期望。

本发明实施例中，解码器模块将后验分布p_θ(y′ⁱ|xⁱ，ε)表述为解码网络，即decoder，由于解码网络decoder相当于原因样本xⁱ和加性噪声N到重构后的结果y′ⁱ的映射模型，所以我们采用与编码网络相同的网络结构。在解码网络decoder中，利用加性噪声N和原因样本xⁱ对结果样本yⁱ进行重构，得到重构后的结果变量y′ⁱ。然后通过计算因果关系中的结果变量yⁱ与解码网络decoder重构后的结果变量来估计重构误差εⁱ，进一步，变分下界的期望可以采用蒙特卡罗方法进行估计。

更进一步地，上述实施例中，判别器模块将简单分布和后验分布p_θ(y′ⁱ|xⁱ，ε)表述为判别网络，判别网络采用两层全连接网络和一层不带有非线性函数处理的输出层作为网络结构，通过判别网络区分加性噪声来自于简单分布还是后验分布p_θ(y′ⁱ|xⁱ，ε)，并使简单分布向后验分布p_θ(y′ⁱ|xⁱ，ε)靠近。

本发明实施例中，在判别器模块中，定义了一个关于和p_θ(ε)的判别器T(xⁱ，yⁱ；ε)，判别器的网络结构由两层全连接网络和一层不带有非线性函数处理的输出层构成。通过判别网络尽可能的区分出加性噪声ε是来自于目前的推断模型还是来自于p_θ(ε)独立分布，同时迫使向分布p_θ(ε)靠近。

更进一步地，上述实施例中，判别网络的目标函数为：

式中，σ(t)＝(1+e^-t)^-t代表Sigmoid函数，T^*(X，Y；ε)代表目标函数的最优值。

本发明实施例中，通过使用编码器模块、解码器模块和判别器模块的交替处理，来优化边缘对数似然函数对应的变分下界使其收敛，并通过变分下界来确定带有中间变量的因果数据之间的因果方向。

本发明实施例另一方面提供一种级联医疗观测数据的因果推断系统，下文描述的该系统可以与上文描述的方法相互对应参照。

请参阅图4，该系统包括：

数据获取模块410，用于获取级联医疗观测数据，并从级联医疗观测数据中抽取第一变量和第二变量；

模型建立模块420，用于以因果关系中的原因、级联结构中每个深度对应的中间变量以及因果关系中的结果为参数，建立改进后的级联非线性加性噪声模型；

函数构造模块430，用于针对改进后的级联非线性加性噪声模型，构造边缘对数似然函数；

函数分解模块440，用于对边缘对数似然函数进行对抗变分解，并利用近似后验分布方法优化得到边缘对数似然函数对应的变分下界；

参数求解模块450，用于以第一变量作为因果关系中的原因，以第二变量作为因果关系中的结果，利用预置的对抗训练模型求解最大化的变分下界，得到第一变分下界值；以第二变量作为因果关系中的原因，以第一变量为作为因果关系中的结果，利用对抗训练模型求解最大化的变分下界，得到第二变分下界值；

方向确定模块460，用于将第一变分下界值与第二变分下界值进行比较，得到比较结果，并根据比较结果，确定级联医疗观测数据的因果方向。

以上可知，本发明实施例提供的级联医疗观测数据的因果推断系统，以因果关系中的原因、级联结构中每个深度对应的中间变量以及因果关系中的结果为参数，建立改进后的级联非线性加性噪声模型，可以更好地匹配具有级联结构的医疗观测数据，提高识别级联医疗数据因果方向的精确度，同时，通过预置的对抗训练模型求解最大化的边缘对数似然函数对应的变分下界，利用对抗策略而非近似公式来绕开KL散度，能够允许加性噪声有更广泛的分布，从而提升模型的推断能力，与现有技术相比，能够很好地识别具有级联结构的医疗观测数据的因果方向，显著提高因果方向识别的精确度。

下面通过具体的实施案例来说明通过本发明实施例公开的级联医疗观测数据的因果推断方法及系统的实际应用结果。

以推断“胰岛素含量-食量-体重”此因果方向为例，在正确的因果关系中，胰岛素含量为初始原因，体重为最终结果，食量为初始原因和最终结果之间的中间变量。

首先，从医疗观测数据中抽取胰岛素含量数据和体重数据。

其次，以胰岛素含量为原因X，体重为结果Y，和随机变量u输入进编码器模块encoder，通过三个带有ReLU非线性函数的全连接层和一个不做非线性处理的输出层网络。

将编码器模块encoder的输出和用户自定义噪声ε输入进判别器，通过两层全连接网络和一层不带有非线性函数处理的输出层网络，使判别器目标函数：得到最优值。

将编码器模块encoder的输出和胰岛素含量X一起输入进解码器模块decoder，通过三个带有ReLU非线性函数的全连接层和一个不做非线性处理的输出层网络。计算出重构误差通过编码器模块、判别器模块、解码器模块交替处理，使得变分下界ELBO收敛，并计算出正向的变分下界值。

再次，改变输入，以体重为原因X，胰岛素含量为结果Y，通过同样的方法，计算出反向的变分下界值。

最后，将正向的变分下界值与反向的变分下界值进行比较，比较出正向的变分下界值大于反向的变分下界值，从而推断出正确的因果方向：“胰岛素含量-食量-体重”。

本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。还需要说明的是，在本说明书中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括该要素的过程、方法、物品或者设备中还存在另外的相同要素。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其他实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。