从时态数据中发现因果关系：概述和新视角

摘要。

时态数据代表了复杂系统的时间顺序观察，一直是工业、医学和金融等许多领域广泛生成的典型数据结构。 分析此类数据对于各种应用程序都非常有价值。 因此，在过去几十年中，人们提出了不同的时态数据分析任务，例如分类、聚类和预测。 其中，因果发现，即从时间数据中学习因果关系，被认为是一项有趣但关键的任务，并引起了广泛的研究关注。 根据时间数据是否经过校准，即，现有的因果发现工作可分为多变量时间序列因果发现和事件序列因果发现两个高度相关的类别。 然而，以往的大多数调查只关注时间序列的因果发现，而忽略了第二类。 在本文中，我们详细说明了这两个类别之间的相关性，并对现有解决方案进行了系统概述。 此外，我们还提供公共数据集、评估指标和时态数据因果发现的新视角。

1. 介绍

记录复杂系统状态变化的时态数据被不同的应用领域广泛收集，例如社交网络、生物信息学、神经科学和金融等。作为最流行的数据结构之一，时态数据由按时间排序的属性序列组成。 由于传感器和计算设备的快速发展，时态数据分析的研究工作近年来不断涌现。 针对不同的任务提出了不同的方法，例如分类（Ismail Fawaz 等人，2019；Ratanamahatana 和 Keogh，2004）、聚类（Aghabozorgi 等人，2015；Liao，2005）、预测（Weigend，2018）、因果发现（Edinburgh 等人，2021；Krakovská 等人，2018）、等。

在这些任务中，识别许多时间成分之间的因果关系的因果发现已成为时态数据分析的一项具有挑战性但又至关重要的任务。 学习到的因果结构可能有助于解释数据生成过程和指导数据分析方法的设计。 根据数据是否经过校准，用于因果发现的时间数据可分为两类，即，多变量时间序列（MTS）和事件序列。 因此，现有的因果发现方法也可以分别分为两类。 在本次调查中，我们的目标是提供深思熟虑的概述并总结时态数据因果发现的前沿。

MTS 数据描述了随时间变化的多个变量的校准状态，是许多领域中的通用时态数据。 从 MTS 中​​发现因果关系可能有利于数据分析模型的可解释性和稳健性。 然而，因果关系的定义并不唯一，导致解决方案不同。 因此，现有研究可分为四类：即基于约束的方法、基于分数的方法、基于函数因果模型（FCM）的方法和基于格兰杰因果关系的方法。 此外，还存在塔肯斯因果关系、微分方程等观点。 在本文中，我们将详细说明每个类别的主要思想和最新进展。

本次调查讨论的另一个任务是从事件序列中发现因果关系，它推断不规则和异步观察的时间序列中的因果关系。 具体来说，它以一系列不同事件作为输入，并输出一个表示不同事件之间因果关系的因果图。 这项任务非常重要，因为大多数现实世界事件不能在固定的时间间隔内出现。 根据MTS任务，我们将相应的方法分为三大类：基于约束的方法、基于分数的方法和基于格兰杰因果关系的方法。 在这三类中，基于格兰杰因果关系的方法，特别是基于格兰杰因果关系的霍克斯过程模型，由于格兰杰因果关系和霍克斯过程之间存在自然匹配而得到了很好的发展。 我们将在本次审查中进一步详细描述这些方法。

表格1。 关于因果发现的现有评论的要点。

Reviews	Multivariate Time-series					Event Sequence	Highlights
	Constrain-based	Score-based	FCM-based	Granger	Deep Learning
(Glymour et al., 2019)	No	No	No	Yes	No	No	An overview for causal discovery methods with practical issues and insightful guidelines
(Guo et al., 2021)	No	No	No	No	No	No	Causal discovery methods dealing with big data (high-dimensional, mixed data) are reviewed
(Vowels et al., 2023)	No	No	No	Yes	No	No	A more extensive coverage of continuous optimization approaches compared to other surveys
(Chen et al., 2022b)	No	No	No	No	No	No	A wider concept of deep learning causal discovery methods is introduced
(Moraffah et al., 2021)	Yes	No	Yes	Yes	No	No	The first survey covers the current progress to analyze time series from a causal perspective
(Shojaie and Fox, 2021)	No	No	No	Yes	Yes	No	Recent advances including network-form and more general notions of Granger causality
(Assaad et al., 2022b)	Yes	Yes	Yes	Yes	No	No	A recent and comprehensive review for causal discovery in time series with comparative evaluations
Ours	Yes	Yes	Yes	Yes	Yes	Yes	A systematic review of causal discovery in both MTS and event sequence, with new perspectives

最近，许多调查（Glymour 等人，2019；Guo 等人，2021；Vowels 等人，2023；Chen 等人，2022b；Moraffah 等人，2021；Shojaie 等人，2021；Assaad 等人，2022b； Kitson 等人，2021；Deng 等人，2022b；Heinze-Deml 等人，2018） 已发表，总结了因果发现的进展。 我们在表1中比较了代表性评论及其亮点。 如图所示，这些调查分为两类。 第一行的研究工作（Glymour 等人，2019；Vowels 等人，2023；Guo 等人，2021；Chen 等人，2022b）从不同角度讨论了一般因果发现问题。 例如，(Glymour 等人, 2019) 简要回顾了计算因果发现方法。 (Vowels 等人, 2023) 关注持续优化方法的快速发展。 为了处理大数据，（郭等人，2021）介绍了基于机器学习的因果推理和因果发现方法。 此外，深度学习因果发现方法在不同的变量范式中进行了综述（陈等人，2022b），从更广泛的角度讨论数据中的因果关系。 在这些论文中，时态数据被视为一种特殊应用，许多数据特定的方法并未包括在内。 第二行的调查重点是时间数据因果发现。 如表1所示，（Edinburgh 等人，2021；Krakovská 等人，2018）对双变量时间序列的因果发现方法进行了综述。 （Moraffah 等人，2021；Shojaie 和 Fox，2021） 最近对时间序列中的因果推理方法进行了综述。 最近的工作（Assaad等人，2022b）讨论并比较评估了时间序列因果发现的现有解决方案。 然而，这些评论忽略了事件序列的因果发现方法。 在本文中，我们不仅对两种时态数据的因果发现方法进行了深入的概述，而且还分析了它们之间的联系和差异。

接下来，我们首先介绍第2节中因果发现问题的背景和初步知识。 MTS 和事件序列因果发现的最新进展分别在 3 和 4 节中详细说明。 之后，我们在第 5 节中概述了时态数据因果发现的应用，并在第 6 节中总结了可用资源。 最后，我们在7节中讨论了最新时态数据因果发现方法的局限性和新视角。 本次调查的整体框架如图1所示。

对于tree=grow'=east，anchor=west，节点选项=draw，thick，font=，align=center，，edge=semithick，分叉边缘，l sep=8mm，s sep=8mm ，文字宽度=2.3cm，叉距=2mm，，[时间因果发现，fill=col1，parent，rotate=90，font=，对于tree=s sep=2.0mm，[MTS

因果发现

(§ 3，表3)，字体=, for tree=child, fill=col3, text width = 4.0cm, [Constraint-Based Approaches
(§ 3.1), text width = 3.6cm, [With Causal Sufficiency,fill = col3,text width = 2.8cm [ oCSE (Sun et al., 2015), PCGCE (Assaad et al., 2022a), PCMCI (Runge et al., 2019b; Runge, 2020), text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] [Without Causal Sufficiency,fill = col3,text width = 2.8cm [ ANLTSM (Chu and Glymour, 2008), tsFCI (Entner and Hoyer, 2010), SVAR-FCI (Malinsky and Spirtes, 2018), LPCMCI (Gerhardus and Runge, 2020), text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] ] [Score-Based Approaches
(§ 3.2), text width = 3.6cm, for tree = child, fill = col3,text width = 3.6 cm [Combinatorial Search,fill = col3,text width = 2.8cm [ Structural EM (Friedman et al., 1998), Greedy Hill-climbling Search (Peña et al., 2005), Structural Constraints (de Campos and Ji, 2011), etc., text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] [Continuous Optimization,fill = col3,text width = 2.8cm [ DYNOTEARS (Pamfil et al., 2020), NTS-NOTEARS (Sun et al., 2021), IDYNO (Gao et al., 2022), text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] ] [FCM-Based Approaches
(§ 3.3), for tree = child, fill = col3,text width = 3.6 cm [Independent Component Analysis,fill = col3,text width = 2.8cm [ VAR-LiNGAM (Hyvärinen et al., 2008, 2010a), MCD (Schaechtle et al., 2013), NCDH (Wu et al., 2022b), text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] [Additive Noise Model,fill = col3,text width = 2.8cm [ TiMINo (Peters et al., 2013), NBCB (Assaad et al., 2021), text width=5.5 cm,draw=colline1,line width=1.2pt,fill=col1,] ] ] [Granger Causality
Based Approaches
(§ 3.4), for tree = child, fill = col3,text width = 3.6cm [HSIC-Lasso-GC (Ren et al., 2020), (R)NN-GC (Montalto et al., 2015; Wang et al., 2018), MPIR (Wu et al., 2020), NGC (Tank et al., 2022), eSRU (Khanna and Tan, 2020), SCGL (Xu et al., 2019), GVAR (Marcinkevics and Vogt, 2021), TCDF (Nauta et al., 2019), CR-VAE (Li et al., 2023), InGRA (Chu et al., 2020), ACD (Löwe et al., 2022), etc., text width=9.4 cm,draw=colline1 ,line width=1.0pt,fill=col1, ] ] [Others
(§ 3.5), for tree = child, fill = col3,text width = 3.6cm [Information-theoretic Statistics (Schreiber, 2000; Runge et al., 2012a; Sun and Bollt, 2014), Differential Equation Based Methods (Voortman et al., 2010; Bellot et al., 2022), Nonlinear State-space Methods (Sugihara et al., 2012), Logic-based Methods (Kleinberg and Mishra, 2009), Hybrid Methods (Li et al., 2016), etc., text width=9.4 cm,draw=colline1 ,line width=1.0pt,fill=col1, ] ] ] [Event Sequence Causal Discovery
(§ 4), font=, for tree=child, fill=col2, text width = 4.0cm [Multivariate Point Process
(§ 4.1), for tree = child, fill = col2, text width=3.6cm, [Basics: Intensity Function, Log-likelihood, text width=9.4 cm,draw=colline1 ,line width=1.0pt,fill=col1, ] ] [Granger Causality
Based Approaches
(§ 4.2), for tree = child, fill = col2,text width=3.6cm [GLM Point Process, text width=2.8cm, [ GLM Model (Kim et al., 2011), text width=5.5cm,draw=colline2,line width=1.2pt,fill=col1] ] [Hawkes Process, text width=2.8cm, [MLE-SGLP (Xu et al., 2016), THP (Cai et al., 2021), $L_0$Hawkes (Idé et al., 2021), HGEM (Yu et al., 2020), NPHC (Achab et al., 2017), GC-nsHP (Chen et al., 2022a), MDLH (Jalaldoust et al., 2022), etc., text width=5.5cm,draw=colline2,line width=1.0pt,fill=col1] ] [Wold Process, text width=2.8cm, [Granger-Busca (de Figueiredo et al., 2018), VI-MWP (Etesami et al., 2021), text width=5.5cm,draw=colline2,line width=1.0pt,fill=col1] ] [Neural Point Process, text width=2.8cm, [CAUSE (Zhang et al., 2020), text width=5.5cm,draw=colline2,line width=1.0pt,fill=col1] ] ] [Others
(§ 4.3), for tree = child, fill = col2,text width=3.6cm [Constraint-Based Approaches, text width=2.8cm, [ MMP-LR/NI (Bhattacharjya et al., 2022), CA (Meek, 2014), text width=5.5cm,draw=colline2,line width=1.2pt,fill=col1] ] [Score-Based Approaches, text width=2.8cm, [PGEM (Bhattacharjya et al., 2018), text width=5.5cm,draw=colline2,line width=1.0pt,fill=col1] ] ] ] [Applications
(§ 5, Table 5),font=, for tree=child, fill=col5, text width = 4.0cm, [Scientific Endeavors, for tree = child, fill = col5,text width = 3.6cm, [ Earth Science, Neuroscience, Bioinformatics, etc., text width=9.4cm,draw=colline3,line width=1.0pt,fill=col1] ] [Industrial
Implementations, for tree = child, fill = col5,text width = 3.6cm, [Anomaly Detection, Root Cause Analysis, Business Intelligence in Online Systems, Video Analysis, Urban Data Analysis, Clinical Data Analysis, etc., text width=9.4cm,draw=colline3,line width=1.0pt,fill=col1] ] ] [Discussions &
New Perspectives
(§ 7) , font=, for tree=child, fill=col4, text width = 4.0cm, [Challenges & Practical Considerations
(§ 7.1), for tree = child, fill = col4, text width=3.6cm, [1) Non-stationarity, 2) Heterogeneity, 3) Unobserved Confounders, 4) Subsampling, 5) Expert Knowledge, text width=9.4 cm,draw=colline4 ,line width=1.0pt,fill=col1, ] ] [New Perspectives
(§ 7.2), for tree = child, fill = col4,text width = 3.6cm [1) Amortized Paradigm, 2) Supervised Paradigm, 3) Causal Representation Learning, text width=9.4cm,draw=colline4,line width=1.0pt,fill=col1] ] ] ]

图1。 从时态数据中发现因果关系的框架。

2. 背景和基本知识

本节首先定义因果发现中的关键概念和假设，然后概述适用于时态数据的三种因果图表示。 最后，将提出从 MTS 和事件序列发现因果关系的问题定义。

表 2. 本次调查中使用的主要符号。

Notation	Description
$d,E$	number of time-series variate, and of event types, respectively
$x_i^t$	the $i$-th time series at time $t$ in multivariate time series
$N_e(t)$	the number of the event $e$ occurrences before time $t$
$⟂⟂,\cancel{}⟂⟂$	independent, and not independent
$V,U$	the set of endogenous variables, and of exogenous variables, respectively
$𝒢$	causal graph
$Pa(x_i)$	the parent nodes of $x_i$

2.1. 因果发现中的关键概念和假设

一些关键概念是从时间数据推断因果关系的基础。 我们在讨论研究工作之前建立这个共同点。 随后，我们提出了结构因果模型、$d$-分离、因果马尔可夫条件、因果可识别性和因果最小性的正式定义，其符号详见表2。

图 2. 基本 DAG 和简单的结构因果模型。

结构因果模型（SCM）。 （Pearl，2009） 中提出的 Pearl 综合因果理论使我们能够利用因果层次 (PCH) 从观察中得出因果结论（Pearl 和 Mackenzie，2018） 。 由此，结构因果模型被定义为因果关系的图形表示，它捕获对一个或多个变量的干预如何影响数据生成机制中其他变量的值。 形式上，SCM可以用4元组表示，其中$V,U$分别表示内生变量和外生变量的集合，$P(U)$是外生变量的分布，$F$表示映射函数的集合。 具体来说，对于 $f_i\in F$，模型 $x_i:=f_i(Pa(x_i),u_i),i=1,\mathrm{\dots },d$ 表示将值 $x_i$ 分配给其结构父项 $Pa(x_i)$ 和外生变量的函数$u_i$。 对于每个 SCM，我们可以通过为每个 $x_i$ 添加一个顶点并从 $Pa(x_i)$ 中的每个父变量（原因）定向边来生成因果图 DAG $G$给子$x_i$（效果）。 SCM与对应的DAG的关系如图2(a)(b)所示。

d-分隔符。 $d$-分离是确定图形模型中两组变量之间不存在因果效应的标准。 如果两组变量之间的每条路径都被阻塞，则称这两组变量是 $d$ 分隔的。 在形式上，一组变量 $𝐒$ d - 如果 $𝐒$ 阻止两个变量之间的所有路径，则将它们分隔开。 对于图 2 (d)(e)(f) 中给定的因果图，两个顶点 $x_1,x_3$ 由顶点集 d 分隔$𝐒$如果$x_2\in 𝐒$。 对于图2(c)（又名v-struct或collider)中的关系，$x_1,x_3$为如果 $x_2$ 且 $x_2$ 的后代都不在集合 $𝐒$ 中，则也用 d 分隔。

$d$-分离是因果发现中的一个基本概念，因为它提供了确定两组变量是否因果相关的标准。 如果两组变量是$d$分开的，那么它们之间不存在直接或间接的因果效应，并且在给定观察到的变量的情况下可以认为它们是独立的。 相反，如果两组变量不是 $d$ 分隔的，那么它们之间可能存在直接或间接的因果效应，在从数据推断因果关系时需要考虑这一点。 因此，$d$-分离是识别图模型中因果关系的重要工具。

因果马尔可夫条件。 在 SCM 的因果图中，每个变量与其直接原因的非影响无关（Pearl，2009）。 换句话说，一个变量在给定其父项（即直接导致它的变量）的情况下有条件地独立于其非影响（即不直接导致它的变量）。 这个条件在因果推理中起着至关重要的作用。 它能够从非实验数据中识别因果效应。 形式上，因果马尔可夫条件意味着联合分布可以根据以下分解进行因式分解：

$$P(𝐱)=\prod _i^{d}P(x_i|Pa(x_i))$$

马尔可夫等价类（MEC）。 如果两个图模型在观测变量之间具有相同的条件独立关系集，则它们属于同一个 MEC，无论图的具体结构如何。 例如，图2(d)(e)(f)中的因果图暗示着相同的d-分离信息$x_1⟂⟂x_3|x_2$并且属于相同的MEC。 MEC 很重要，因为它使我们能够识别从非实验数据推断因果关系所需的最小条件集。

因果可识别性。 如果可以在不做出任何不可检验的假设或调用观察变量之外的附加信息的情况下进行估计，则因果效应是可识别的。 这意味着从观察的角度来看，同一 MEC 中的所有因果图都代表等效的因果结构。 一般来说，因果可识别性要求因果图是非循环的，并且治疗和结果变量之间的所有后门路径都被阻止。 如果满足这些条件，则可以使用$do$-演算或其他因果推理技术来识别因果效应。 因此，因果发现的前提是因果关系是可识别的。

因果极小性。 考虑 DAG $𝒢$ 和概率分布 $P$,$P$ 被认为满足相对于 $𝒢$ 的因果极小性，如果 $P$ 相对于 $𝒢$ 是马尔可夫的，但不是 $𝒢$ 的任何真子图。 它表明所有变量对于准确表示因果关系都是必要且充分的，同时排除任何对因果机制无贡献的变量。 当且仅当在给定剩余父节点的情况下，不存在条件独立于其任何父节点的节点时，相对于因果图的分布是最小的。 换句话说，所有的父母都是“主动”的。

基于上述概念，我们引入了三个假设：因果充分性、忠实性和时间优先性，它们是因果发现不可检验的基础。

因果充分性。 如果观察到所有变量的所有共同原因，则一组变量是因果充分的（Spirtes等人，2000）。 这一假设表明，SCM中的因果图能够反映真实数据的生成过程，并且不存在隐藏的混杂因素。 在因果充分性的假设下，大多数因果发现算法都假设因果结构可以描述为 DAG。

忠诚。 忠实断言观测数据中存在的所有 $P$ 条件独立关系都由因果模型 $𝒢$ 蕴含，反之，因果模型隐含的所有条件独立关系也成立在观察到的数据中。 请注意，忠实意味着因果最小性。 如果 $P$ 相对于 $𝒢$ 是忠实且马尔可夫的，则满足因果极小性。

图 3. 忠实性和因果最小性的例子。

直觉上，忠诚并不容易理解。 我们尝试用一个例子（Peters等人，2017）来阐明它。 如图3所示，我们假设${𝒢}_1$的生成过程为线性高斯SCM：

$X$ $:=N_X$

$Y$ $:=aX+N_Y$

$Z$ $:=bY+cX+N_Z$

噪声变量$N_X\sim 𝒩(0,{\sigma }_x^2)$、$N_Y\sim 𝒩(0,{\sigma }_y^2)$和$N_Z\sim 𝒩(0,{\sigma }_z^2)$是共同独立的。 让我们考虑一个特殊情况，即$a\cdot b+c=0$。 在此设置中，变量 $X$ 和 $Z$ 是独立的。 $Y\to Z$ 的方向将反转，因果模型 ${𝒢}_1$ 退化为 ${𝒢}_2$。 根据定义，${𝒢}_1$和${𝒢}_2$满足因果极小性。 但在这种特殊情况下，忠实性被违反了，即即${𝒢}_2$不是${𝒢}_1$的真子图。 因此，该线性高斯模型的概率对于${𝒢}_1$并不忠实。 尽管 ${𝒢}_2$ 是 $\mathscr{H}$ 的真子图，但分布不满足因果极小性，因为相对于 $\mathscr{H}$ 的概率不是马尔可夫分布。

虽然忠实度在实践中是无法测试的，但它对于从数据中得出有效的因果推论至关重要，因为它确保模型正确地表示数据生成机制。 如果违反这个假设，因果关系是不确定的，这对于因果发现方法来说是一场灾难（Spirtes等人，2000）。

时间优先。 对于两个变量，时间优先意味着原因必须先于其结果发生。 它是从时间数据中发现因果关系的基础假设，并在因果过程中创建了不对称的时间关系。 当两个变量存在因果关系时，时间优先级有助于我们建立因果关系的方向。 然而，如果时间序列的采样频率较高，则与时间序列相关的事件之间的时间差可能难以区分。 在这种情况下，在不同时间发生的两个事件在观察时间序列中可以被视为瞬时发生，从而导致在不同时刻发生的原因和结果之间存在同时的因果关系。

2.2. 时态数据的因果结构

图 4. 离散时间观察的因果图。

对于时间数据，因果关系可以通过时间优先级（Eichler，2012）直观地定义，表明原因先于其影响。 它揭示了时间上不对称的因果关系，当已知两个变量存在因果关系时，可用于确定因果关系。 基于时间优先级，因果结构存在三种图形表示，即、全时因果图、窗口因果图和总结因果图。

如图4（a）所示，全时因果图代表了动态系统的完整图。 对于$d$变量​​时间序列$𝐱$，每个时间点$t$的测量是一个向量$(x_1^t,\mathrm{\dots },x_d^t)$。 全时因果图中的顶点由每个时间点 $t$ 处的组件 $x_1,\mathrm{\dots },x_d$ 集以及特定于滞后的有向链接（例如 $x_i^{t-k}\to x_j^t$)组成。 然而，由于每个系列在每个时间点的单一观察，通常很难发现全时因果图。

为了解决这个问题，提出了窗口因果图。 它假设了一种时间同构因果结构，即观测向量 $𝐱$ 的动态受 的支配，其中函数 $f$ 根据过去的 和噪声 ${𝐮}^t$ 决定接下来的观测。 如图4（b）所示，窗口因果图在时间窗口中表示，其大小相当于全时因果图中的最大滞后。

如图4（c）所示，每个时间序列组件被折叠成一个节点以形成汇总因果图。 汇总图表示时间序列之间的因果关系，不涉及时间滞后（Peters等人，2013）。 在许多应用中，无需精确了解时刻之间的相互作用，就足以对时间变量之间的关系进行建模。

对于从时间数据中发现因果关系，大多数工作的目标是找到总结因果图。 然而，总结因果图并不总是对应于 SCM，这意味着它们无法实现与底层时间解析 SCM 一致的干预预测（Janzing 等人，2018；Rubenstein 等人，2017）。

2.3. 问题定义

图 5. 问题描述：(a) 从 MTS 数据发现因果关系，(b) 从事件序列数据发现因果关系。

如图5所示，时态数据的因果发现可以分为两个问题，即，MTS的因果发现和事件序列的因果发现。 接下来，我们分别正式定义它们。

MTS 的因果发现。 考虑一个带有 $d$ 变量的时间序列：${\{{𝐱}^t\}}_{t\in {\mathbb{Z}}^{+}}={\{{(x_1^t{x}_2^t\mathrm{\dots }{x}_d^t)}^{\top }\}}_{t\in {\mathbb{Z}}^{+}}$。 假设变量之间的因果关系由以下结构方程模型给出：

$$x_i^t:=f_i(Pa(x_i^t),u_i^t),i=1,\mathrm{\dots },d,$$

其中对于任何 $i\in \{1,\mathrm{\dots },d\}$ 在时间实例 $t$ 中，$Pa(x_i^t)$ 是 $x_i^t$ 的直接父母集，可以是过去的，也可以是同时实例。 $u_i^t$表示独立噪声，可以表示测量噪声或驱动噪声（Peters等人，2022）而不失一般性。 MTS 的因果发现旨在找到两种输出中的任何一种，即，总结因果图或窗口因果图。 对于汇总因果图，输出是总结因果结构的邻接矩阵$A\in {\mathbb{Z}}^{d\times d}$，矩阵$A$的第$(i,j)$项为$1$如果$x_i$过去的观察进入$x_j^t$的结构方程，否则进入$0$的结构方程。 如果$A_{ij}=1$，我们就说“$x_i$导致$x_j$”。 对于最大时间滞后$K$的窗口因果图，输出矩阵$W$和$A^k(k\in \{1,\mathrm{\dots },K\})$分别对应于片内和片间边缘。 例如，$W_{ij}=1$ 表示瞬时依赖$x_i^t\to x_j^t$，而$A_{ij}^k=1$ 表示$k>0$ 的滞后依赖$x_i^{t-k}\to x_j^t$。

从事件序列中发现因果关系。 对于事件序列：$\{(t_1,e_1),(t_2,e_2),\mathrm{\dots }\}$、$t_i,i=1,2,\mathrm{\dots }$表示事件发生的时间，而$e_i,i=1,2,\mathrm{\dots }$表示对应的事件类型。 我们的目标是发现不同事件类型之间的因果关系。 一般来说，我们可以构造一个因果图$G=(g_i),i=1,2,\mathrm{\dots },n$，其中每个节点代表一种事件序列。 我们的使命是发现因果图中的边缘。 例如，如果从节点$g_j$到节点$g_i$有一条有向边，我们说事件类型$g_j$是事件类型$g_i$。

3. 从多元时间序列发现因果关系

在本节中，我们回顾多元时间序列数据的因果发现方法，包括基于约束的方法、基于评分的方法、基于函数因果模型的方法、格兰杰因果关系等。 代表性算法结合其特点总结于表3中。

表3。 按类别排列的时间序列数据的因果发现算法的特征。

Section	Method	Causal Graph	Nonlinear	Instantaneous effects	Hidden confounders	Sufficiency Asm.	Markov Asm.	Faithfulness Asm.	Minimality Asm.
Constraint-based	oCSE (2015) (Sun et al., 2015)	Summary	Yes	No	No	Yes	Yes	Yes
	PCGCE (2022) (Assaad et al., 2022a)	Extended	Yes	Yes	No	Yes	Yes	Yes
	PCMCI (2019) (Runge et al., 2019b)	Window	Yes	No	No	Yes	Yes	Yes
	PCMCI+ (2020) (Runge, 2020)	Window	Yes	Yes	No	Yes	Yes	Yes
	ANLTSM (2008) (Chu and Glymour, 2008)	Window	Yes	Yes	Yes	No	Yes	Yes
	tsFCI (2010) (Entner and Hoyer, 2010)	Window	Yes	No	Yes	No	Yes	Yes
	SVAR-FCI (2018) (Malinsky and Spirtes, 2018)	Window	No	Yes	Yes	No	Yes	Yes
	FCIGCE (2022) (Assaad et al., 2022a)	Extended	Yes	Yes	Yes	No	Yes	Yes
	LPCMCI (2020) (Gerhardus and Runge, 2020)	Window	Yes	Yes	Yes	No	Yes	Yes
Score-based	DYNOTEARS (2020) (Pamfil et al., 2020)	Window	No	Yes	No	Yes	Yes	No	No
	NTS-NOTEARS (2021) (Sun et al., 2021)	Window	Yes	Yes	No	Yes	Yes	No	No
	IDYNO (2022) (Gao et al., 2022)	Window	Yes	Yes	No	Yes	Yes	No	No
FCM-Based	VAR-LiNGAM (2008) (Hyvärinen et al., 2008)	Window	No	Yes	No	Yes	Yes	No	Yes
	NCDH (2022) (Wu et al., 2022b)	Summary	Yes	No	No	Yes	Yes	No	Yes
	TiMINo (2013) (Peters et al., 2013)	Summary	Yes	Yes	No	Yes	Yes	No	Yes
	NBCB (2021) (Assaad et al., 2021)	Summary	Yes	Yes	No	Yes	Yes	Yes	Yes
Granger Causality	HSIC-Lasso-GC (2020) (Ren et al., 2020)	Summary	Yes	No	No	No	No	No	No
	(R)NN-GC (2015,2018) (Montalto et al., 2015; Wang et al., 2018)	Summary	Yes	Yes	No	No	No	No	No
	MPIR (2019) (Wu et al., 2020)	Summary	Yes	No	No	No	No	No	No
	NGC (2022) (Tank et al., 2022)	Summary	Yes	No	No	No	No	No	No
	eSRU (2020) (Khanna and Tan, 2020)	Summary	Yes	No	No	No	No	No	No
	SCGL (2019) (Xu et al., 2019)	Summary	Yes	No	No	No	No	No	No
	GVAR (2021) (Marcinkevics and Vogt, 2021)	Summary	Yes	No	No	No	No	No	No
	TCDF (2019) (Nauta et al., 2019)	Window	Yes	Yes	Yes	No	No	No	No
	CR-VAE (2023) (Li et al., 2023)	Summary	Yes	Yes	No	No	No	No	No
	InGRA (2020) (Chu et al., 2020)	Summary	Yes	No	No	No	No	No	No
	ACD (2022) (Löwe et al., 2022)	Summary	Yes	No	Yes	No	No	No	No
Others	DBCL (2010) (Voortman et al., 2010)	Summary	Yes	Yes	Yes	No	Yes	Yes
	NGM (2022) (Bellot et al., 2022)	Summary	Yes	Yes	No	No	No	No	No
	CCM (2012) (Sugihara et al., 2012)	Summary	Yes	No	No	No	No	No	No
	PCTL(c) (2009,2011) (Kleinberg and Mishra, 2009; Kleinberg, 2011)	Summary	Yes	No	No	No	No	No	No

3.1. 基于约束的方法

作为一系列因果发现算法，基于约束的方法依赖于条件独立性的统计测试，并且易于理解和广泛使用。 我们首先给出基于约束的方法的主要思想，包括一般步骤和因果假设。 详细的方法将分为有因果充分性假设和无因果充分性假设的方法，并分别介绍。

总体步骤是：首先，基于条件独立性构建变量之间的骨架。 其次，根据规则中的定向标准对骨架进行定向。 目标是构建代表真实因果图 MEC 的完整部分有向无环图 (CPDAG)。 这些从观察中得出 MEC 的方法的核心是因果假设。 这些方法通常基于因果马尔可夫性质和忠实度的假设，有些还假设因果充分性（没有未观察到的混杂因素）。 在本节中，我们首先回顾主要算法及其对假设因果充分性的时间序列数据的扩展，然后介绍在因果充分性假设无法保证的情况下的方法。

3.1.1. 因果充分性方法

在这一部分中，我们回顾因果充分性的方法。 为了揭示这些方法的原理，我们首先简要介绍非时间环境中的方法。 然后，基于两种类型的扩展（转移熵和瞬时条件独立性检验），回顾了几种流行的基于约束的时间序列方法，这些方法源于非时间数据的方法。

对于从非时间数据中提取因果关系，Sprites-Glymour-Scheines (SGS) 算法(Spirtes 等人, 1990) 是最早的基于约束的方法之一，被证明是一致的在假设因果充分性的独立同分布（i.i.d）观察下。 然而，它面临着耗尽所有节点之间的独立性测试的问题。 非常大的搜索问题使其不适合实际。 Peter-Clark (PC) 算法(Spirtes 等人, 2000) 也假设因果充分性，其引入是为了减少不必要的条件独立性测试和搜索过程。 给定$d$非时间变量，PC算法的详细过程定义如下3步： (1)首先，算法从完整的无向图$G$开始。 (2) 其次，当$n=0,1,2,\mathrm{\dots },d-2$时，算法分别检索是否存在以其他$n$变量为条件的变量对$i$和$j$。 。 如果满足，则删除$i$和$j$之间的无向边，并将条件变量更新为分离集。 它继续到修剪后的骨架。 (3)最后确定获得CPDAG的碰撞体(V结构)，并根据其他规则确定剩余的无向边。

尽管 SGS 和 PC 等方法是在非时间设置中设计的，但基于约束的时间序列数据方法通常是从它们扩展而来的。 我们将回顾最近四种流行的基于约束的方法，这些方法也假设时间序列数据的因果充分性。 其中，两个扩展(Sun 等人, 2015; Assaad 等人, 2022a)基于传递熵的因果概念，另外两个(Runge 等人, 2019b; Runge, 2020）其中通过瞬时条件独立性测试扩展到时间序列。

基于传递熵的时间序列扩展。 传统的基于约束的方法可以基于传递熵的概念扩展到时间序列的场景。 传递熵是一种无模型的时间因果关系度量，其定义和变体将在3.5.1小节中详细介绍。 在这里，我们将转移熵度量视为现成的部分，并从基于约束的方法论的角度回顾两种代表性方法。

最优因果熵(oCSE)原理（Sun等人，2015）被提出来指导从MTS数据中计算和数据高效的因果发现算法。 它基于因果熵的理论概念，是传递熵的推广，用于测量许多变量的成对关系到网络关系 oCSE 方法采用的过程与 PC 中的过程略有不同：而不是尽可能限制其条件的大小设置后，它从一开始就以构成所有可用节点的过去的所有潜在原因为条件。 该算法总结为算法1，它由因果节点的聚合发现和非因果节点的逐步去除组成。 具体来说，给定节点$j$，联合执行两个过程来推断其直接因果邻居：（1）首先，它发现$j$的超集$Pa(x_j)$'基于因果熵最大化的直接因果邻居聚合。 (2)其次，根据因果熵准则修剪掉非直接因果邻居，例如，如果$\mathrm{CE}(x_i^t\to x_j^{t+1}|Pa(x_j^t)\backslash \{x_i^t\})=0$，则从$Pa(x_j)$中删除$i$。 这是一种计算效率高、样本效率高的算法。 然而，它假设隐藏的动态遵循平稳的一阶马尔可夫过程，因为因果熵仅模拟时滞等于 1 的因果关系。 最近，提出了PCGCE (Assaad 等人, 2022a)，基于PC算法和贪婪因果熵来提取时间序列数据的扩展汇总因果图，是因果熵的一个变体。

Input: Multivariate time series $𝐱$ with $d$ dimensions, a significant threshold $\alpha$

Output: The summary causal graph $G$

1 Initialize an empty graph $G$ with $d$ nodes $V$

2 for $j\in \{1,\mathrm{\dots },d\}$ do

3  # Aggregative Discovery of Causal Nodes

4  $z=\mathrm{\infty }$

5  while $z>0$ and card($Pa(x_j)$) do

6   for $x_i\in V\backslash Pa(x_j)$ do

7    Compute the p-value ($z_p$) corresponding to the test $\mathrm{CE}(x_i^t\to x_j^{t+1}|Pa(x_j^t))>0$

8   if $z_p>\alpha$ then add edge $x_i\to x_j$ to $G$

9  # Progressive Removal of Non-Causal Nodes

10  for $x_i\in Pa(x_j)$ do

11   Compute $z$ corresponding to the test $\mathrm{CE}(x_i^t\to x_j^{t+1}|Pa(x_j^t)\backslash \{x_i^t\})=0$

  if $z>\alpha$ then remove edge $x_i\to x_j$ from $G$

Algorithm 1 oCSE

通过瞬时条件独立性测试扩展到时间序列。 PCMCI算法（Runge等人，2019b）利用PC算法的变体，灵活地结合线性或非线性条件独立性测试，并从时间序列数据中提取因果关系。 该算法的目标是发现窗口因果图。 与PC算法不同，PCMCI首先构建一个部分连通图，其中所有节点对$(x_i^{t-k},x_j^t)$如果$k>0$则定向为$x_i^{t-k}\to x_j^t$。 此初始化还迎合时间优先级。 该算法由两个阶段组成：(1) 与 PC 中所做的一样，PCMCI 基于条件独立性删除所有不必要的边缘。 此外，它还基于时间一致性的假设来删除同源边缘。 (2)利用瞬时条件独立(MCI)来处理自相关，这可能会导致虚假相关。 这里，MCI 是一种测量，在测试 $X_i^{t-k}\cancel{}⟂⟂X_j^t|Pa(X_j^t)\mathrm{\backslash }\{X_i^{t-k}\},X_i^{t-k}$ 时以 $x_j^t$ 和 $x_i^{t-k}$ 的父母为条件。 它还提供了从 $x_i^{t-k}$ 到 $x_j^t$ 因果强度的可解释概念。 PCMCI 已被证明是一致的，并且可以灵活地与任何类型的条件独立性测试（线性或非线性）结合，例如部分相关性和互信息。 近年来，也有大量针对非参数检验的机器学习方法，解决了广泛的独立性和依赖性类型（Zhang等人，2011；Runge，2018）。

PCMCI⁺ 算法（Runge，2020） 扩展了 PCMCI，包括瞬时因果关系的发现。 PCMCI⁺ 算法的核心是两个偏离原始 PC 算法的基本思想：首先，它分别对滞后条件集和同时期条件集进行边缘去除过程。 其次，它利用 MCI 来校准自相关下的 CI 测试，这与 PCMCI 中的类似。 作者在（Runge，2020）中还详细介绍了自相关的祸与福。

3.1.2. 没有因果充分性的方法

本部分将回顾没有因果充分性的基于约束的方法。 首先，我们简单介绍一下非时间数据的快速因果推理（FCI）算法（Spirtes等人，2000）。 那么，针对MTS数据的方法包括两类：（1）通过时间序列模型进行快速因果推理，它是从FCI扩展而来的。 (2) 通过瞬时条件独立性测试的方法。

FCI算法是PC算法的推广，可以在存在潜在混杂因素的情况下使用，并被证明是渐近正确的。 它利用对观测数据的独立性测试来提取有关观测变量之间祖先关系的（部分）信息，因此 FCI 算法的目标是推断适当的 PAG。 FCI 算法首先构建一个由无向边组成的完整图，类似于 PC 算法。 然后进行迭代条件独立性测试以去除边缘。 因此，FCI 算法首先在使用 Sepset 和 Possible-Dsep 集进行调节时删除独立的边。 对于剩余的无向边，递归地应用十个方向规则。 详细的FCI算法，包括理论分析，证明了算法的合理性和完整性，可以在（Zhang，2008）中找到。

通过时间序列模型进行快速因果推断。 基于约束的方法称为加性非线性时间序列模型(ANLTSM)（Chu and Glymour，2008）是在假设隐藏混杂因素的影响是线性且同时发生的情况下提出的。 为了避免非参数条件独立性检验的维数灾难，ANLTSM 利用加性回归模型，可以指定如下：

$$x_j^t=\sum _{1\le i\le d,i\ne j}{a}_{j,i}{x}_i^t+\sum _{1\le i\le d,1\le l\le \tau }{f}_{j,i,l}(x_i^{t-l})+\sum _{r=1}^h{b}_{j,r}{u}_r^t+e^t.$$

这里，$a_{j,i}$和$b_{j,r}$是常数值，$f_{j,i,l}(\cdot )$表示平滑单变量函数。 多维高斯白噪声形式的未观测效应可分为两类：$e^t$反映观测变量的潜在直接原因，${(u_r^t)}_{1\le r\le h}$表示潜在的共同原因。 潜在的共同原因同时影响观察到的变量。 对于$x_i^t$和$x_j^t$,$u_r^t$足以被表述为潜在的共同原因，当且仅当存在$1\le r\le h$使得$b_{j,r}{b}_{i,r}\ne 0$。 基于上述的加性回归模型，利用FCI算法来识别滞后和瞬时因果关系。 为了检测瞬时关系，首先通过估计条件期望$𝔼(x_i^t|x_j^t\cup S)$，在给定集合$S$的情况下测试$x_i^t$和$x_j^t$之间的条件独立性，然后使用 F 检验或 BIC 分数等统计检验来检查 $x_i^t$ 和 $x_j^t$ 之间的预测关系的显着性，其中预测变量的不显着性意味着之间的条件独立性$x_i^t$ 和 $x_j^t$。 滞后因果关系以类似的方式识别。 其余的边根据规则定向。 如果数据生成满足加性非线性时间序列模型，则该方法被证明是一致的。 然而，ANLTSM 方法将同期相互作用限制为线性，并将潜在混杂因素限制为线性和同期。

FCI 对时间序列数据的另一个扩展是 tsFCI （Entner 和 Hoyer，2010） 算法，其中 FCI 算法直接通过时间窗口应用。 具体来说，假设观察到的时间序列数据来自处于平衡状态的系统，原始时间序列数据通过大小为$\tau$的滑动窗口转换为一组随机向量样本。 然后将变换后的向量的每个分量视为单独的随机变量，直接应用原始FCI算法。 由于从标准 FCI 导出的信息量非常有限，因此进一步将时间优先级和时间不变性作为背景知识纳入其中，以便在定向阶段做出更多推断。 然而，tsFCI 忽略了选择变量和同期因果关系。 最近，提出了一种名为 SVAR-FCI （Malinsky 和 ​​Spirtes，2018） 的基于约束的方法，该方法允许数据生成过程中的瞬时影响和任意潜在混杂。 与 tsFCI 类似，它也使用时间不变性来推断额外的边缘去除。

通过瞬时条件独立测试的方法。 研究发现，由于条件独立性检验的效应量较小（Gerhardus 和 Runge，2020），原始 FCI 算法及其时间变体在自相关时间序列情况下召回率较低。 一些研究人员旨在在存在未观察到的混杂变量的情况下扩展 PCMCI，以解决上述问题。 在（Gerhardus and Runge，2020）中，提出了Latent PCMCI(LPCMCI)算法。 LPCMCI算法的核心有两个思想：第一，基于因果发现中效应大小的分析，它使用变量的父代作为默认条件，并且在条件集中不测试非祖先，这不仅避免了夸大的错误积极但也减少了要测试的集。 其次，它引入了中间标记和 LPCMCI-PAG 的概念，作为明确的因果解释，以促进边缘的早期定位。 LPCMCI算法被证明是阶无关的、健全的和完整的。

3.2. 基于分数的方法

另一类因果发现方法是基于评分函数的。 首先将介绍基于评分的方法的主要思想，包括（动态）贝叶斯网络、基于评分的方法与基于约束的对应方法相比的特征、模型评分和模型搜索。 然后，我们将分别回顾 MTS 的组合搜索方法和连续优化方法。

3.2.1. 基于分数的方法的基础知识

基于分数的方法的动机是，编码错误（条件）独立性的图结构也会导致模型拟合不佳。 在基于分数的方法中，因果结构附加到贝叶斯网络（BN）或动态贝叶斯网络（DBN）的概念（Dean和Kanazawa， 1989；Murphy，2002） 处理时间数据。 鉴于此，基于评分的方法可以生成多个模型并对其进行概率评分，然后输出最可能的模型。 这与基于约束的方法形成鲜明对比，后者派生并输出单个模型，而没有量化其正确的可能性。 通过应用拟合优度测量而不是条件独立性检验，在基于评分的方法中，忠实性假设被稀释。 因此，从观察中学习 BN 或 DBN 的问题可以表述为：给定一组实例，找到最匹配它们的网络，即优化目标函数。 它由两个元素组成：模型评分和模型搜索。

模型评分。 常见的目标函数分为两类：贝叶斯分数，重点关注拟合优度并允许结合先验知识；信息论分数，明确考虑模型复杂性，除了拟合优度之外，还旨在避免过度拟合-of-fit (Kitson 等人, 2021)。 贝叶斯评分函数家族包括贝叶斯狄利克雷等效（BDe）评分（Heckerman等人，1995）、K2评分（Kayaalp和Cooper，2013）等。 最广泛使用的信息论评分包括贝叶斯信息准则 (BIC) （Neath 和 Cavanaugh，2012 年） 和 Akaike 信息准则 (AIC) （Burnham 和 Anderson，2004 年） t1>.

模型搜索/优化。 基于分数的方法将搜索因果结构$G$的问题转化为使用上述分数函数$S$的优化程序。因此最终目标表述为(Peters等人, 2017)：

$$\widehat{G}={\mathrm{argmin}}_{G\mathrm{over}𝐱}S(D,G),$$

其中$D$表示变量$𝐱$的经验数据。 传统上，这是一个组合图搜索问题，并且解决方案通常不是最优的，因为众所周知，找到全局最优网络是 NP 困难的（Chickering，1995）。 由于图的搜索空间很大，贪婪等价搜索（GES）（Chickering，2002）等一系列工作涉及局部启发式。 然而，它们仍然遭受维数诅咒和次优问题的困扰。 最近，表征非循环约束的代数结果被用于结构学习，将组合问题转化为连续优化问题（Zheng等人，2018，2020），可以重新表述为：

其中 $𝐀$ 表示邻接矩阵，$h$ 是用于在推断结构中强制非循环性的函数。 原始无环约束函数实现为(Zheng 等人, 2018)中的$h(𝐀)=\mathrm{tr}(e^{𝐀\odot 𝐀})-d$。 它依靠增强拉格朗日方法（ALM）（Yurkiewicz，1985）来解决连续约束优化问题。 各种工作进一步采用神经网络中的连续约束公式来提取非线性因果关系（郑等人，2020；于等人，2019；高等人，2021）。

在时间序列的背景下，基于分数的方法的最终目标是学习 DBN 的结构。 DBN是一种以时间序列为变量的概率网络，它可以分解为先验网络和过渡网络。 先前网络提供给定时间戳中变量之间的依赖性，而过渡网络提供随时间变化的依赖性。 因此，DBN 代表同一框架中的同期效应和延时效应。 基于对时间序列的扩展，我们回顾了基于分数的方法，遵循从组合搜索到连续约束优化的类似范式。

3.2.2. 组合搜索方法

为了有效地从MTS数据中进行基于评分函数的组合搜索，研究人员开发了多种方法，包括结构期望最大化（Friedman等人，1998）、交叉验证（Peña等人， 2005），以及评分函数的分解（de Campos and Ji，2011）。

在(Friedman等人, 1998)中，作者首先利用了Structural Expectation-Maximization (Structural EM)算法(Friedman, 1997, 1998) ，最初是推断BN的标准算法，从纵向数据中学习DBN。 结构 EM 算法将结构和参数修改与单个 EM 过程相结合，可以找到由评分函数定义的局部最优值。

在（Peña等人，2005）中，$K$折叠交叉验证（CV）被用作学习DBN的计算上可行的评分标准。 给定观察数据$D$，该数据被随机分割成大约相等大小的$K$折叠$D^1,\mathrm{\dots },D^K$，模型$G$的CV值公式为$\frac{1}{T}{\sum }_{k=1}^K\log p(D^k|G,{\widehat{\theta }}^k)$。 并使用贪婪爬山搜索来估计$E[\log p(D_{T+1}|G,\widehat{\theta })]$。 该过程从空图开始，并通过应用可用的最高得分单边附加或删除来逐渐更新它。 实验表明，对于大范围的样本量，基于交叉验证的评分方法比基于 BDe 的 BIC 的模型具有更好的泛化能力。

基于可分解的评分函数，论文（de Campos and Ji，2011）利用结构约束将DBN中的结构学习问题转化为相应的增广BN，并提出了一种分支和-bound 算法保证全局最优性。 最优目标的分解形式可以形式化为：

$$(G^{0\ast },G^{\prime \ast })={\mathrm{argmax}}_{G^0,G^{\prime }}(S_{D_0}(G^0)+S_{D_{1:T}}(G^{\prime }))=({\mathrm{argmax}}_{G^0}{S}_{D_0}(G^0)+{\mathrm{argmax}}_{G^{\prime }}{S}_{D_{1:T}}(G^{\prime })),$$

其中$G^0$和$G^{\prime }$分别对应于先前网络和过渡网络。 结构约束作为减少搜索空间的一种方式，指定可以或不可以包含弧的位置。 由于分支定界特性，算法可以在当前最佳解和全局最优解的上限处停止。 受益于分支定界算法和结构约束，所提出的方法被证明能够处理比以前更大的数据集。

3.2.3. 持续优化方法

由于NOTEARS （郑等人，2018）最近的贡献，DAG的基于分数的学习可以重新表述为连续约束优化问题，这启发了各种工作（郑等人） , 2020; Yu 等人, 2021; Ng 等人, 2022b, 2020) 结构学习。 该方法的核心是表示为约束函数的非循环性的代数表征，进一步利用该代数特征来最小化最小二乘损失，同时强制执行非循环性。 在时间序列的背景下，一些作品也采用了这种连续约束公式来支持结构学习和因果发现（Pamfil等人，2020；Sun等人，2021；Hsieh等人，2021；Gao等人，2022 ）。

DYNOTEARS，在（Pamfil等人，2020）中引入，通过连续优化方法从时间序列数据中捕获线性关系。 它采用以下标准 SVAR 方式对数据进行建模：

$${𝐱}^t={𝐱}^t𝐖+{𝐱}^{t-1}{𝐀}^1+\mathrm{\dots }+{𝐱}^{t-p}{𝐀}^p+{𝐮}^t,$$

其中$p$是SVAR模型的阶数，$𝐮$是中心误差变量的向量。 为了保证 SVAR 模型的可识别性，误差项 ${𝐞}^t$ 被假设为非高斯或标准高斯，即 ${𝐮}^t\sim 𝒩(0,I)$，因为可识别性被证明在两种情况下都成立 （Hyvärinen 等人，2010b；Peters 等人，2017）。 $𝐖$和$𝐀$是加权邻接矩阵，分别对应片内边缘（同期关系）和片间边缘（时滞关系）。 SEM还可以采用紧凑的形式：${𝐗}^t={𝐗}^t𝐖+{𝐗}^{(t-p):(t-1)}𝐀+𝐔$。 结构学习的过程围绕最小化受非循环约束的最小二乘损失，这给出了以下优化问题：

${\min }_{𝐖,𝐀}f(𝐖,𝐀)\mathrm{s}.\mathrm{t}.$ $𝐖\mathrm{is}\mathrm{acyclic},$

$\mathrm{where}f(𝐖,𝐀)={\displaystyle \frac{1}{2n}}||{𝐗}^t-{𝐗}^t𝐖-{𝐗}^{(t-p):(t-1)}𝐀$ $|{|}_F^2+{\lambda }_{𝐖}||𝐖|{|}_1+{\lambda }_{𝐀}||𝐀|{|}_1.$

为了避开非循环约束下求解优化问题的关键难点，DYNOTEARS 遵循（郑等人，2018）中的工作，其中迹指数函数$h(𝐖)=\mathrm{tr}(e^{𝐌\odot 𝐌})-d$被利用为无环性的等效表述。 连续约束优化问题通过增强拉格朗日方法转化为以下形式的无约束问题：

$${\min }_{𝐖,𝐀}F(𝐖,𝐀),\mathrm{where}F(𝐖,𝐀)=f(𝐖,𝐀)+\frac{\rho }{2}h{(𝐖)}^2+\alpha h(𝐖).$$

针对上述平滑增强目标的优化，分别提出了两种求解方法。 第一种方法是使用标准求解器，例如 L-BFGS-B (Zhu 等人, 1997)。 另一种方法是类似于 (Hyvärinen 等人, 2010b) 中的两阶段过程，我们可以将方程重写为 $𝐳=\mathrm{𝐳𝐖}+𝐔$ 并得出 $𝐖$ 通过对错误项 $𝐳$ 使用静态注释。

NTS-NOTEARS (Sun 等人, 2021) 是采用连续约束公式的最新进展。 与线性自回归模型 DYNOTEARS 相比，NTS-NOTEARS 能够提取变量之间的线性和非线性关系。 它通过利用一维卷积神经网络 (CNN) 来实现这一目标，该网络利用输入数据中的顺序拓扑，因此非常适合时态数据的神经函数近似模型。 $d$ CNN，其中每个第一层都是具有 $m$ 内核的一维卷积层，联合训练，其中 $j$ CNN 预测目标的期望给定先前和同期输入变量的特定时间 $t$ 的变量 $x_t^j$。 每个 CNN 都可以被视为目标变量的马尔可夫毯子。 DBN 中子变量对其父变量的依赖关系如下：

$$𝔼[x_j^t|Pa(x_j^t)]={\mathrm{CNN}}_j(\{{𝐱}^{t-k}:1\le k\le K\},{𝐱}_{-j}^t),$$

其中父级 $Pa(x_j^t)$ 来自经过训练的 CNN，${𝐱}_{-j}^t$ 表示时间步 $t$ 处除 $x_j$ 之外的所有变量。 根据NOTEARS-MLP(Zheng 等人, 2020)（NOTEARS的非线性和基于神经网络的扩展），DBN中边的依赖强度通过以下方式估计：

$$W_{ij}^k={\Vert {\varphi }_{i,j}^k\Vert }_L^2\mathrm{for}k=1,\mathrm{\dots },K+1.$$

具体地，当估计的依赖强度大于阈值权重$W_{ij}^k>W_{thres}^k$时，$x_i^{t-k}$属于父集合$Pa(x_j^t)$。 优化过程遵循与 DYNOTEARS 类似的方式。 还值得注意的是，NTS-NOTEARS 显示了变量依赖性的先验知识，这些知识可以转换为额外的优化约束并合并到 L-BFGS-B 求解器中。

为了同时处理观察数据和介入数据，最近提出了一种名为 IDYNO (Gao 等人, 2022) 的算法。 它首先通过神经网络引入非线性目标来模拟复杂的动态，然后修改客观和通用的解决方法来处理干预目标的不同分布。

我们可以发现，使用持续优化并避免可能因果结构的显式组合遍历，是基于分数的结构学习的一种强大方法。 过去几年也见证了这种方法的大量应用和扩展。 然而，(Kaiser and Sipos, 2022; Reisach 等人, 2021; Ng 等人, 2022a)进一步讨论了一些边界和限制，包括数据规模的影响和增强模型的收敛条件。拉格朗日方法。 我们建议您在进一步开发和应用该方法系列时考虑这些问题。

3.3. 基于 FCM 的方法

上述两种方法要么面临 MEC 的不可分离性，要么需要大样本来确认因果忠实性。 因果发现也可以基于功能因果模型（FCM）（Pearl等人，2000）进行，在2.1中也称为SCM，通过以下方式描述因果系统：一组方程。 近年来，基于 FCM 的时态和非时态数据方法的激增。 在本小节中，我们首先介绍基于 FCM 的方法的主要思想，包括函数因果模型和使用噪声来定向因果关系。 然后，将分别评述基于 FCM 的两个方法系列，即使用独立分量分析和加性噪声模型的方法。

在 FCM 中，每个变量都通过方程根据其直接原因和一些附加噪声进行解释。 例如，函数$x_j=f_j(x_i,u_j)$用一些额外的噪音$u_j$解释了因果关系$x_i\to x_j$。 基于 FCM 的因果发现方法的一个基本思想是，统计噪声可以成为一种有价值的洞察来源，它迎合了最近的发现（Climenhaga 等人，2021）的挑战性正统观念认为噪音应被视为滋扰。 具体来说，可以借助噪声来识别和估计因果关系。

3.3.1. 使用独立成分分析的方法

在这一部分中，我们首先通过回顾非时间环境下的原始算法（Shimizu等人，2006）来介绍这一系列方法的基本思想。 然后详细介绍MTS数据的方法(Hyvärinen等人, 2008, 2010a; Schaechtle等人, 2013; Wu等人, 2022b)。

LiNGAM (Shimizu 等人, 2006)是一种典型的非时间环境下基于 FCM 的因果发现算法，具有以下假设：（1）线性数据生成过程，（2）非高斯扰动，(3) 没有未观察到的混杂因素。 在LiNGAM模型中，观测值之间的关系可以表示为$𝐱=\mathrm{𝐁𝐱}+𝐮$，其中$𝐱,𝐁,𝐮$分别表示变量向量、因果图邻接矩阵和噪声向量。 该等式可以重写为$𝐱=\mathrm{𝐀𝐮}$，其中$𝐀={(𝐈-𝐁)}^{-1}$。 对于方程，可以使用独立成分分析（ICA）方法（Stone，2004）来估计$𝐀$和因果关系$𝐁$可以推导出来。 沿着这个思路，DirectLiNGAM (Shimizu 等人, 2011) 进一步利用回归模型来确保原始模型以受控的步骤数收敛到正确的解决方案。 LiNGAM 对时间序列的扩展如下。

作为 LiNGAM 的时间扩展，VAR-LiNGAM (Hyvärinen 等人, 2008, 2010a) 利用非高斯性质来估计结构自回归 (SVAR) 模型。 SVAR 模型反映瞬时和时滞因果效应，是实证经济学中分析动态现象最流行的工具（Moneta 等人，2013）。 在VAR-LiNGAM中，时间序列的表示是SVAR和SEM的组合，定义为：

(SVAR) $${𝐱}^t=\sum _{k=0}^{\tau }{𝐁}^k{𝐱}^{t-k}+{𝐮}^t$$

其中 ${𝐁}^k$ 是变量 $𝐱$ 之间具有时间滞后 $k$ 的因果效应的 $n\times n$ 矩阵。 ${𝐮}^t$ 是对外部影响或“干扰”进行建模的随机过程，假定它们是独立的、时间上不相关且非高斯的。 为了估计上述模型，结合了自回归（AR）模型（时间滞后$k>0$)的经典最小二乘估计，其形式化为：

(VAR) $${𝐱}^t=\sum _{k=1}^{\tau }{𝐌}^k{𝐱}^{t-k}+{𝐧}^t$$

基于SVAR和VAR形式化，VAR-LiNGAM的基本思想是我们可以以经典的最小二乘方式一致且高效地估计VAR模型的${𝐌}^k$。 通过LiNGAM分析我们可以推导出瞬时因果效应的估计。 至于延时效应，可以通过重参数化推导出来。 接下来的方法详细定义如下，分四个步骤： (1)首先，拟合回归并用${\widehat{𝐌}}^k$表示AR矩阵的最小二乘估计。 (2) 其次，计算残差，即${\widehat{𝐧}}^t={𝐱}^t-{\sum }_{k=1}^{\tau }{\widehat{𝐌}}^k{𝐱}^{t-k}$。 (3) 第三，根据方程${\widehat{𝐧}}^t={𝐁}^0{\widehat{𝐧}}^t+{𝐞}^t$进行LiNGAM分析（Shimizu等人，2006），得到瞬时因果效应的估计${\widehat{𝐁}}^0$。 (4) 最后，计算时滞因果效应${\widehat{𝐁}}^k(k>0)$的估计值${\widehat{𝐁}}^k=(𝐈-{\widehat{𝐁}}^0){\widehat{𝐌}}^k$。 如果自回归部分的阶数设置为零，即$\tau =0$，则VAR-LiNGAM模型将退化为LiNGAM模型。 这种方法在实证经济学中的深入应用可以在（Moneta等人，2013）中找到。

将VAR-LiNGAM扩展到时变情况下因果模型的识别和估计（黄等人，2015），其中进一步利用高斯过程回归来自动建模因果模型如何变化时间。 在(Lanne等人, 2017)中，初始VAR-LiNGAM被推广到推断图可以包含环的情况。 并且所提出的模型在理论上被证明是可识别的。 （Schaechtle等人，2013）中提出了另一种基于LiNGAM的算法，称为多维因果发现（MCD）。 MCD可以通过集成数据分解和投影来有效地发现多维设置（例如时间序列数据）中的因果依赖关系。

为了摆脱线性(Hyvärinen 等人, 2008, 2010a) 或加性假设(Peters 等人, 2013) 的约束，提出了一种基于 FCM 的算法，名为“非线性因果发现”最近在（Wu等人，2022b）中提出了通过HM-NICA(NCDH)来从时间序列中提取一般非线性关系。 该算法的核心是利用非线性 ICA 算法来测量非线性关系。 假设观察结果是由相互独立的潜在组件生成的：

$$𝐱=𝐟(𝐒)\mathrm{where}𝐟={(f_1,f_2,\mathrm{\dots },f_d)}^T\mathrm{and}𝐒={(S_1,S_2,\mathrm{\dots },S_d)}^T.$$

与线性 ICA 类似，$𝐒$ 包含彼此独立的分量，非线性 ICA 的目标是从 $𝐱$ 恢复 $𝐒$。 NCDH 首先利用非线性 ICA 结合 HMM （Hälvä 和 Hyvärinen，2020） 来分离潜在噪声。 作为对 ICA 排列不确定性的补救措施，进行了一系列独立性测试以确定观测变量与分离噪声之间的对应关系。 最后采用递归搜索算法来提取因果关系。

3.3.2. 使用加性噪声模型的方法

事实上，有许多非线性因果关系违反了 LiNGAM 系列方法的假设。 尽管最近在一般非线性条件下提取因果关系取得了进展（例如 NCDH），但它们的使用受到限制。 另一类基于 FCM 的方法基于具有非线性函数的加性噪声模型 (ANM)，适用于更一般的设置。 在这一部分中，首先给出使用ANM的方法的主要思想。 然后我们将介绍MTS数据的详细方法。

(Hoyer等人, 2008)证明，如果因果极小条件成立，则可以在具有非线性函数的ANM中识别真实的因果结构。 在ANM中，如果$x_i\to x_j$，我们有$x_j=f(x_i)+u_j$，并且原因$x_i$和加性噪声$u_j$是独立的。 如果噪声$u$服从非高斯分布且$f(\cdot )$是线性函数。 在双变量情况$x_i\to x_j$中，我们可以在因果和反因果方向上拟合回归模型，通过残差检验独立性可以推断真实方向。 对于多变量情况，可以采用成对策略（Mooij等人，2009）。 (Peters等人, 2014)讨论了该算法的正确性。

(Peters等人, 2013)提出了独立噪声时间序列模型(TiMINo)，这是一种基于ANM的时间序列因果发现方法。 它输入时间序列数据并输出汇总时间图或保持未定状态，这避免了在模型指定错误或数据不足时导致错误的因果结论。 它利用与非时间和多变量设置中类似的方法（Mooij等人，2009）。 具体来说，它试图拟合时间序列的结构方程模型，其可以表述如下：

$$x_j^t=f_j(Pa{(x_j^{\tau })}^{t-\tau },\mathrm{\dots },Pa{(x_j^1)}^{t-1},Pa{(x_j^0)}^t,u_j^t),$$

其中误差项 $u_j^t$ 与变量索引 $j$ 和时间索引 $t$ 共同独立。拟合方法有多种选择$f$，例如线性模型、广义加性模型和高斯过程回归模型。 为了推断加性噪声模型中的因果关系，可以利用互相关和 HSIC 等独立性测试（Gretton 等人，2007）。

这些函数因果模型（例如 VAR-LiNGAM 和 TiMINo）存在一些缺点。 结果表明，这些方法随着节点数量的增加而无法很好地扩展（Glymour等人，2019），并且在没有大样本的情况下，这些方法的性能并不乐观（Malinsky和Danks，2018） ）。 为了克服这些缺点，（Assaad 等人，2021）中提出了基于噪声/基于约束(NBCB)的方法，其中基于约束的方法进一步基于时间序列数据的原始加性噪声模型。 具体来说，每个时间序列的潜在原因都是通过类似于TiMINo中的加性噪声​​模型来检测的。 使用时间因果熵来修剪不必要的因果关系，时间因果熵是因果熵（Sun等人，2015）的扩展，用于测量两个时间序列之间的（条件）依赖性。

3.4. 基于格兰杰因果关系的方法

格兰杰因果关系是在许多实际应用中分析时间序列数据的流行工具。 存在许多基于格兰杰因果关系开发的因果发现方法。 在本小节中，我们首先介绍格兰杰因果关系的定义。 在深入研究详细方法之前，将给出并比较 MTS 的两类格兰杰因果关系模型（无模型和基于模型）。 由于基于模型的方法在更一般的情况下具有优越性，本部分的其余部分将重点关注基于模型的方法的两个最新进展：（1）基于核的方法(3.4.3) ，（2）基于神经网络的方法(3.4.4)。

3.4.1. 格兰杰因果关系的基础知识

格兰杰因果分析最早于（Granger，1969）提出，是一种基于可预测性确定因果关系的强大方法。 如果$x_i$的过去值提供了关于$x_j$的未来值的唯一的、统计上显着的信息，则时间序列$x_i$格兰杰原因$x_j$。 根据这个命题，如果 $x_i$ 被定义为 $x_j$ 的“因果”，如果

其中表示给定所有相关信息的历史记录的$x_j^t$的最佳预测。这里表示从中排除的信息。上述定义看似笼统，没有具体的建模假设，但为了更好的表示能力和方便推理，格兰杰因果关系也有多种形式的定义，基于不同的模型规范和统计工具，例如自回归模型（格兰杰原文）论文（Granger，1969）)等等。 如果观察到所有相关变量并且不存在瞬时联系，则格兰杰因果关系相当于底层 DAG 中的因果关系(Peters 等人, 2013, 2017)。

3.4.2. MTS 的早期方法

早期识别格兰杰因果关系的方法仅限于二元设置。 具体来说，双变量设置中格兰杰因果分析的一个有据可查的（Lütkepohl，1982）问题是，如果不调整所有相关协变量，因果发现可能会产生误导。 一方面，有必要考虑更多变量，以防止识别错误的格兰杰因果关系（Shojaie and Fox，2021）。 另一方面，MTS广泛存在于各个领域。 在MTS中推断格兰杰因果关系，在一些文献中也称为图形格兰杰因果关系或网络格兰杰因果关系，已成为研究热点。 MTS 的各种图形 Granger 因果分析模型可以分为两类，即无模型方法和基于模型的方法。

无模型方法。 多元格兰杰因果关系的无模型方法的主流是基于可预测性，并且需要估计条件概率密度函数（CPDF）（Bai等人，2010）。 在（Diks and Wolski，2016）中，提供了CPDF的估计，并将双变量Diks-Panchenko非参数因果关系检验扩展到多变量情况。 通过在边际概率密度函数中引入条件变量，基于copula的格兰杰因果关系模型（Hu和Liang，2014；Kim等人，2020）也可以扩展到多元情况。 此外，转移熵和定向信息等无模型度量（Amblard 和 Michel，2011） 能够检测非线性依赖性。 这些无模型估计器的定义和一些属性将在3.5.1中详细介绍。 无模型方法可以很好地处理非线性格兰杰因果关系。 然而，这些估计器存在高方差，需要大量数据才能进行可靠估计，并且当变量数量增加时，还会遭受维数灾难。 因此，在非线性和高维的复杂现实场景中，无模型方法的使用在某种程度上受到限制。

基于模型的方法。 与无模型方法相比，基于模型的方法计算效率高，因此更适合推断高维条件下的格兰杰因果关系。 绝大多数格兰杰因果模型都采用基于模型的推理方法，其中测量的时间序列通过合适的参数化数据生成模型进行建模。 推断出的参数最终揭示了格兰杰因果关系的真实拓扑。 沿着这条线的早期方法通常在线性时间序列动态假设下使用流行的向量自回归（VAR）模型。 对于$d$变量​​时间序列$𝐱$，VAR模型定义为：

$${𝐱}^t=\sum _{k=1}^{\tau }{A}^k{𝐱}^{t-k}+{𝐮}^t,$$

其中 $A^k$ 是 $d\times d$ 矩阵，指定滞后 $k$ 如何影响序列的未来演变，${𝐮}^t$ 表示零均值噪声。 在 VAR 模型中，作为双变量情况 （Granger，1969） 的直接扩展，时间序列 $i$ 不会导致时间序列 $j$ 如果并且仅当对于所有时间滞后$k$,$A^k$ 的分量$(j,i)$ 都等于0。 因此，格兰杰因果分析减少了确定 $A^k$ 中的哪些条目在所有滞后上为零。 还有大量的研究工作（Arnold等人，2007；Lozano等人，2009a；Shojaie和Michailidis，2010；Basu等人，2015）通过Lasso惩罚及其变体降低计算复杂度高维时间序列中的格兰杰因果分析，也称为套索格兰杰因果关系(Lasso-GC)。 对于这些方法，格兰杰因果序列选择问题一般可以基于最小二乘损失表述如下：

$${\min }_{A^1,\mathrm{\dots },A^{\tau }\in {ℝ}^{d\times d}}\sum _{t=\tau +1}^T{\Vert {𝐱}^t-\sum _{k=1}^{\tau }{A}^k{𝐱}^{t-k}\Vert }_2^2+\lambda R(𝐀),$$

其中 $R(\cdot )$ 是稀疏性诱导正则化器，具有多种实现，如表 4 所示。 不同的惩罚项在 $A^1,\mathrm{\dots },A^{\tau }$ 中引入不同的稀疏模式，从而在格兰杰因果序列选择中引入不同的启发式和约束。 除了 Lasso-GC 之外，另一个值得一提的基于多变量环境下 VAR 模型的工作是条件 Granger 因果指数（CGCI）（Geweke，1982）。 对于变量$X,Y$和条件变量$Z$，通过比较简化模型和完整模型${\mathrm{CGCI}}_{X\to Y|Z}=\ln \frac{\mathrm{var}({ϵ}_{Y|Z})}{\mathrm{var}({ϵ}_{Y|XZ})}$的残差，可以区分多元系统中的直接因果关系和间接因果关系基于CGCI制作。 沿着这一思路，mBTS-CGCI 在（Siggiridou 和 Kugiumtzis，2016） 中提出，基于改进的向后时间选择（mBTS）来限制 VAR 模型的阶数，因此可以更好地应用于高维场景。

表 4. 常见的稀疏性惩罚项，描述为（Nicholson等人，2017；Marcinkevics and Vogt，2021）

Model Structure	Penalty Function
Basic Lasso	${\Vert 𝐀\Vert }_1$
Elastic net	$\alpha {\Vert 𝐀\Vert }_1+(1-\alpha ){\Vert 𝐀\Vert }_2^2,\alpha \in (0,1)$
Lag group Lasso	${\sum }_{k=1}^{\tau }{\Vert {𝐀}^k\Vert }_F$
Component-wise Lasso	${\sum }_{p=1}^d{\sum }_{k=1}^{\tau }{\Vert {({𝐀}^{k:\tau })}_p\Vert }_2$
Element-wise Lasso	${\sum }_{p=1}^d{\sum }_{q=1}^d{\sum }_{k=1}^{\tau }{\Vert {({𝐀}^{k:\tau })}_{p,q}\Vert }_2$
Lag-weighted Lasso	${\sum }_{k=1}^{\tau }{k}^{\alpha }{\Vert {𝐀}^k\Vert }_1,\alpha \in (0,1)$

尽管与无模型方法相比，基于模型的方法可以有效处理高维时间序列，但这些方法的根本问题是模型指定错误。 特别是，基于普通 VAR 模型的多元格兰杰因果关系概念假设时间序列遵循线性动态，而现实世界应用中的许多相互作用本质上是非线性的。 最近，出现了许多与非线性因果关系兼容的基于模型的方法，它们可以分为两类：基于内核的方法和基于神经网络的方法>。 格兰杰因果关系的产生，其基本脉络和发展方向已从经典文献中进行了回顾和展望。 在本小节的以下部分中，由于它们能够在复杂的现实场景中发挥作用，我们将详细介绍基于模型的方法在非线性和高维设置中的最新进展，特别是来自神经网络的新观点。

3.4.3. 基于内核的最新进展

为了在基于模型的方法中提取非线性因果关系，建立非线性参数模型是一种常见的策略。 一系列工作将 Granger 因果关系扩展到核方法（Ancona 等人，2004；Marinazzo 等人，2008b，a；Sindhwani 等人，2013；Ren 等人，2020）。 在(Ancona等人, 2004)中，格兰杰因果关系通过径向基函数扩展到二元非线性情况。 进一步，基于再生核希尔伯特空间（RKHS）理论，提出了格兰杰因果分析模型（Marinazzo等人，2008b）。 关键思想是将数据嵌入希尔伯特空间并搜索该空间中的非线性关系。 然后将该方法推广到(Marinazzo 等人, 2008a) 中的多变量情况。 在（Sindhwani等人，2013）中，提出了核方法的矩阵值扩展，强加于向量值RKHS字典。 该算法针对高维非线性多元回归，自然可以导致图形格兰杰因果关系的非线性推广。 最近，提出了一种基于Hilbert-Schmidt独立准则Lasso Granger因果关系的算法(HSIC-Lasso-GC)（Ren等人，2020）。

3.4.4. 基于神经网络的最新进展

神经网络能够表示变量之间的非线性、复杂和非加性相互作用。 本部分将回顾基于神经网络的格兰杰因果方法的最新进展，包括非均匀嵌入（Montalto等人，2015；Wang等人，2018）、信息正则化（ Wu 等人, 2020)，逐组件神经网络建模 (Tank 等人, 2017, 2022; Khanna 和 Tan, 2020)，低秩逼近 (Xu 等人, 2019)、自解释网络(Marcinkovics and Vogt, 2021)、注意力机制(Nauta 等人, 2019; Schwab 等人, 2019) 、循环变分自动编码器（Li等人，2023）和归纳建模（Chu等人，2020；Löwe等人，2022）。 此外，如图6所示，现有的基于神经网络的格兰杰因果关系方法可以分为四类：基于参数的（Tank等人，2022；Khanna和Tan，2020）、基于注意力的(Nauta 等人, 2019; Chu 等人, 2020)、基于自我解释的(Marcinkovics and Vogt, 2021)、以及关系编码——基于（Löwe等人，2022）。

图 6. 四种基于神经网络的格兰杰因果发现方法根据均值推导因果权重。

具有非均匀嵌入的 DL 扩展。 NN-GC (Montalto 等人, 2015) 中提出了一种称为非均匀嵌入 (NUE) 的特征选择过程，用于识别 MLP 中的显着 Granger 原因模型。 通过贪婪地添加预测时间序列的滞后分量作为输入，MLP 被迭代更新。 如果在过程终止时添加了至少一个滞后分量，则预测时间序列被称为目标时间序列的显着格兰杰原因。 在 RNN-GC (Wang 等人, 2018) 中，通过用门控 RNN 模型替换 MLP 来扩展 NUE，但是，由于该技术需要训练并比较许多候选模型，在高维设置中代价高昂。

具有信息正则化的 DL 扩展。 为了提取非线性动力学，引入了最小预测信息正则化(MPIR)（吴等人，2020）的方法。 它利用预测变量的可学习损坏，并最大限度地减少相互信息正则化风险，它将格兰杰因果关系范式与深度学习模型的优点结合起来。 在 MPIR 中，作者指出，将神经网络与格兰杰因果关系结合起来的简单方法有两个主要缺点：不稳定和低效率。 解决方案是鼓励每个 ${𝐱}_i^{t-K:t-1}$ 向 $x_j^t$ 提供尽可能少的信息，同时通过学习到的损坏保持良好的预测，用简单的方式替换预测 $x_j^t$ 的方式一次丢失一个${𝐱}_i^{t-K:t-1}$。 风险定义如下：

$$R_{𝐗,x_j}[f_{\theta },𝐧]=E_{{𝐗}^{t-1},x_j^t,𝐮}[{(x_j^t-f_{\theta }({\widetilde{𝐗}}_{(𝐧)}^{t-K:t-1}))}^2]+\lambda \cdot \sum _{p=1}^{d}I({\tilde{X}}_{i(n)}^{t-K:t-1};X_i^{t-K:t-1}),$$

其中 ${\widetilde{𝐗}}_{(𝐧)}^{t-K:t-1}:={𝐗}^{t-K:t-1}+𝐧\odot 𝐞$ （或其按元素表示，${\tilde{X}}_{i(n)}^{t-K:t-1}:=X_i^{t-K:t-1}+n_i\cdot e_i,i=1,2,\mathrm{\dots },d$)是具有可学习噪声幅度 $𝐧$ 和 $u_j\sim N(\mathrm{𝟎},𝐈)$ 的噪声损坏输入。 而$W_{pq}=I({\tilde{X}}_{i(n^{\ast })}^{t-K:t-1};X_i^{t-K:t-1})$是$R_{𝐗,x_j}[f_{\theta },𝐧]$最小化时的最小预测信息，包含因果信息，衡量变量$i$对预测变量$j$ 具体来说，$W_{ij}=0$ 如果$x_i⟂⟂x_j$。 此外，由于估计大维度的互信息项效率较低，因此导出了上限作为替代优化目标。 该框架只需要单独训练 $d$ 模型，而不是训练许多候选模型并遭受不稳定和低效率的困扰。

带有组件式神经网络建模的深度学习扩展。 另一种基于神经网络的测量非线性格兰杰因果关系的方法是分量建模。 (Tank等人, 2017)提出了一个组件化的框架，它可以被视为线性VAR模型的推广。 具体来说，各个变量的生成过程可以写成如下：

$${𝐱}_j^t:=g_j({𝐱}_1^{1:(t-1)},\mathrm{\dots },{𝐱}_i^{1:(t-1)},\mathrm{\dots },{𝐱}_d^{1:(t-1)})+u_j^t,\mathrm{for}1\le j\le d,$$

其中 $g_j(\cdot )$ 是基于正则化神经网络实现的连续函数，指定 $𝐱$ 的过去值如何确定变量 ${𝐱}_j$ 的未来值。 在这种情况下，当且仅当 对 不变时，时间序列 对时间序列 ${𝐱}_j$ (${𝐱}_i\nrightarrow {𝐱}_j$) 格兰杰非因果关系，其定义如下："当且仅当 $g_j(\cdot )$ 对 ${𝐱}_i^{1:(t-1)}$ 不变时，时间序列 ${𝐱}_i$ 格兰杰非因果关系：

$$g_j({𝐱}_1^{1:(t-1)},\mathrm{\dots },{𝐱}_i^{1:(t-1)},\mathrm{\dots },{𝐱}_d^{1:(t-1)})=g_j({𝐱}_1^{1:(t-1)},\mathrm{\dots },{𝐱}_{i^{\prime }}^{1:(t-1)},\mathrm{\dots },{𝐱}_d^{1:(t-1)}),$$

对于所有 $({𝐱}_1^{1:(t-1)},\mathrm{\dots },{𝐱}_d^{1:(t-1)})$ 和所有 ${𝐱}_i^{1:(t-1)}\ne {𝐱}_{i^{\prime }}^{1:(t-1)}$。 我们将分别介绍基于该框架的两种方法(Tank 等人, 2022; Khanna and Tan, 2020)。

神经格兰杰因果关系(NGC)在（Tank等人，2022）中提出，使用具有稀疏输入层权重的结构化MLP和LSTM来推断非线性格兰杰因果关系，称为分别为组件级 MLP (cMLP) 和组件级 LSTM (cLSTM)。 在 cMLP 中，每个非线性输出 $g_j$ 均使用单独的 MLP 进行建模，以便轻松理清输入到输出的影响。 第一层的输入矩阵提供格兰杰因果关系的惩罚选择信息。 具体来说，在$g_j(\cdot )$的第一层

$$h_1^t=\sigma (\sum _{k=1}^{\tau }{W}_1^k{𝐱}^{t-k}+b_1),$$

如果权重矩阵$W_1^k$的第$i$列在所有时间滞后$k$中都包含零，则时间序列$i$不会格兰杰-原因系列$j$。与VAR类型方法类似，Granger因果序列通过以下编码选择（Tank等人，2017）过程来选择：

$${\min }_{𝐖}\sum _{t=\tau }^T(x_j^t-g_j(x_{(t-1):(t-\tau )}))+\lambda \sum _{i=1}^{d}R({(W_1)}_{:i}),$$

其中稀疏诱导惩罚$R(\cdot )$是通过组套索惩罚实现的，它提取因果关系而不需要精确的滞后规范。 至于 cLSTM，它回避了滞后选择问题，格兰杰因果信息也可以在普通 LSTM 模型中轻松解释。 输入矩阵与MLP中的略有不同，定义为$W^1={({(W^f)}^{\top },{(W^{in})}^{\top },{(W^o)}^{\top },{(W^c)}^{\top })}^{\top }$，控制过去的时间序列如何影响遗忘门、输入门、输出门和单元更新。 可以根据 $W^1$ 列的组套索惩罚来选择格兰杰因果序列。 最后，为了优化 cMLP 或 cLSTM 中的非凸优化目标，使用近端梯度下降（Parikh 等人，2014），这会导致输入矩阵中精确的零。 优化过程中的这一性质满足了解释框架中格兰杰非因果性的要求。 为了推断格兰杰因果关系的网络拓扑，需要用每个变量作为响应来训练 $d$ 模型。

（Khanna and Tan，2020）中提出了另一种样本高效的架构经济——SRU(eSRU)。 它利用统计循环单位 (SRU)（Oliva 等人，2017） 对观察到的时间序列数据进行建模。 这里，SRU 是一种特殊类型的 RNN，专为具有时滞和非线性依赖性的 MTS 设计，因此也适合提取非线性 Granger 因果关系的网络拓扑。 具体来说，由于非门控架构，它较少受到梯度消失和爆炸问题的影响，并且能够通过维护多时间尺度的汇总统计来对多元时间序列之间的短期和长期时间依赖性进行建模。 与 cLSTM 等基于模型的方法类似，格兰杰因果关系的度量可以从 SRU 的输入层权重参数中导出。 然而，由于因果推理问题中普遍存在的数据稀缺问题，原有的框架存在过拟合的问题。 此外，还实施了两项修改作为 eSRU 中过度拟合的补救措施。

具有低秩近似的 DL 扩展。 （Xu等人，2019）中提出了可扩展因果图学习(SCGL)框架。 作者首先将数据非线性解构为两种类型(即单变量级别非线性和多元级别非线性），并分别建模。 SCGL 的关键思想是，当变量 $d$ 的大小相当大时，学习邻接矩阵 $A\in {ℝ}^{d\times d}$ 的完整大小将是不可扩展的。 实际上，隐藏空间中变量的关系是低秩的（Zorzi and Chiuso，2017；Chiuso and Pillonetto，2012）。 因此，很自然地通过 $k$ 等级分解来近似 $A$，其中 。低秩近似降低了因果发现中的噪声影响，并为下游时间序列分析提供了可解释性（黄等人，2020c）。

具有自解释网络的深度学习扩展。 为了更好的可解释性，提出了广义向量自回归(GVAR)模型（Marcinkovics and Vogt，2021）。 它基于自解释神经网络的扩展（Alvarez-Melis 和 Jaakkola，2018）。 自解释神经网络本质上是由受限属性驱动的可解释模型，并遵循以下形式：

$$f(𝐱)=g(\theta {(𝐱)}_{1}h{(𝐱)}_1,\mathrm{\dots },\theta {(𝐱)}_{k}h{(𝐱)}_k),$$

其中$g(\cdot )$和$𝐡(𝐱)$分别表示链接函数和可解释的基本概念。 结合格兰杰因果推理中经常指定的向量自回归模型，GVAR 模型由下式给出

$${𝐱}^t=\sum _{l=1}^{\tau }{\mathrm{\Psi }}_{{\theta }_l}({𝐱}^{t-l}){𝐱}^{t-l}+{𝐮}^t,$$

其中${\mathrm{\Psi }}_{{\theta }_l}:{ℝ}^d\to {ℝ}^{d\times d}$是由${\theta }_l$参数化的神经网络，其输出是影响力强度对应的矩阵。 具体地，影响$x_i^{t-l}\to x_j^t$的强度是通过${\mathrm{\Psi }}_{{\theta }_l}({𝐱}^{t-l})$的分量$(j,i)$来测量的。 损失函数由三项组成：MSE损失、稀疏诱导正则化（可以从4中选择）和平滑惩罚，定义如下：

$$\frac{1}{T-\tau }\sum _{t=\tau +1}^T{\Vert {𝐱}^t-{\widehat{𝐱}}^t\Vert }_2^2+\frac{\lambda }{T-\tau }\sum _{t=\tau +1}^{T}R({\mathrm{\Psi }}_t)+\frac{\gamma }{T-\tau -1}\sum _{t=\tau +1}^{T-1}{\Vert {\mathrm{\Psi }}_{t+1}-{\mathrm{\Psi }}_t\Vert }_2^2,$$

这里${\{{𝐱}^t\}}_{t=1}^T$是观测到的d变量时间序列，而${\widehat{𝐱}}^t$是GVAR模型做出的一步预测。 现在可以通过 ${\mathrm{\Psi }}_{{\widehat{\theta }}_k}({𝐱}^t)$ 导出每个时间点 $t$ 的解释矩阵，还可以评估格兰杰因果效应的迹象及其随时间的变化。 此外，基于时间反演格兰杰因果关系启发式的GVAR过程（Winkler等人，2016），期望关系在时间反演数据上翻转，以提高稳定性推断的结构。 与上述方法（例如 cMLP、cLSTM、eSRU 和 MPIR）相比，另一个关键区别是这些方法需要训练 $d$ 神经网络，而 GVAR 需要训练 $2\tau$ 网络。

带有注意力机制的深度学习扩展。 (Nauta 等人, 2019) 中引入了时间因果发现框架 (TCDF)，该框架利用基于注意力的扩张 CNN。 该框架由$d$独立的基于注意力的CNN组成，具有相同的架构但不同的目标变量$X_j$。 对于每个目标变量，提出了一个神经网络来导出预测、注意力分数和核权重。 直观上，在预测 $X_j$ 时对 $X_i$ 的关注得分较高，表明前者包含对后者的预测信息。 另外还提​​供了基于排列的程序，用于评估变量重要性并识别重要的因果关系。 TCDF 可以发现自因果关系以及因果关系之间的时间延迟。 此外，通过假设双向因果关系不可能是瞬时的，它还可以以相等的延迟检测隐藏的混杂因素的存在。

此外，IMV-LSTM（Guo等人，2019，2018）提出了一种具有混合注意力的可解释多变量LSTM来提取变量重要性知识。 它被广泛用作多元时间序列中因果发现的基线。 然而，关于注意力及其解释的话题在某种程度上仍然是一个有争议且没有定论的话题（Jain和Wallace，2019；Wiegreffe和Pinter，2019；Grimsley等人，2020）。 特别是在 Granger 因果解释的背景下，天真训练的软注意机制被认为（Sundararajan 等人，2017；Schwab 等人，2019；Chu 等人，2020）不提供任何收益激励准确的归因。 在(Schwab等人, 2019)中，基于平均格兰杰因果误差的度量引入了格兰杰因果注意力权重。 添加 $i$ 时误差的减少可以计算为：$\mathrm{\Delta }{\epsilon }_{X,i}={\epsilon }_{X\mathrm{\backslash }\{i\}}-{\epsilon }_X$，给定辅助预测误差 ${\epsilon }_X,{\epsilon }_{X\mathrm{\backslash }\{i\}}$，无论有没有来自 $i$ 那么格兰杰因果注意力因子可以计算为：${\omega }_i(X)=\frac{\mathrm{\Delta }{\epsilon }_{X,i}}{{\sum }_{j=1}^d\mathrm{\Delta }{\epsilon }_{X,j}}$。 注意因子 ${\omega }_p(X)$ 能够捕获格兰杰因果关系，如果第 $p$ 个时间序列对于目标序列是格兰杰非因果关系，则该因果关系为零。

带有循环变分自动编码器的深度学习扩展。 最近，提出了因果循环变分自动编码器(CR-VAE)（Li等人，2023），其中生成模型将格兰杰因果学习纳入数据生成过程。 通过防止在解码之前对未来信息进行编码，CR-VAE的编码器遵循格兰杰因果关系原理。具体来说，给定时间滞后$\tau$，CR-VAE模型可以写为：

$${\widehat{𝐱}}^{t-\tau :t}=D_{\theta }({𝐱}^{t-\tau :t-1},E_{\psi }({𝐱}^{t-2\tau -1:t-\tau -1}))+{ϵ}^t,$$

其中$E_{\psi },D_{\theta }$代表编码器和解码器。 与经典循环 VAE 的另一个区别是 CR-VAE 利用多头解码器，其中第 $i$ 头被设计用于生成 ${𝐱}_i$。 此外，还利用误差补偿模块来捕获瞬时效应。 CR-VAE 不仅能够提取因果关系，还能够受益于学习到的因果矩阵，以透明的方式进行数据生成过程。

具有归纳建模的深度学习扩展。 归纳建模方法的问题与上述方法略有不同，上述方法收集来自大量个体的 MTS 数据，这些数据具有不同的因果机制，但具有共同的结构。 目标是在具有异质结构的样本上训练模型，以发现每个个体的格兰杰因果关系。 这里回顾了两种归纳建模方法。

(Chu 等人, 2020)提出了归纳格兰杰因果模型(InGRA)，结合Granger因果注意力(Schwab 等人, 2019) 和原型学习。 因为现实世界中经常存在这样的场景，即从具有共同点的异质个体收集大量多元时间序列数据。 InGRA 不是为每个人训练一个或一组模型，而是为可能具有不同格兰杰因果结构的个体训练一个全局模型，避免样本效率低下和过度拟合问题。 首先，利用格兰杰因果注意机制来量化变量对预测的贡献。 由于格兰杰因果注意力不够鲁棒，无法从单个个体的有限数据重建格兰杰因果拓扑，InGRA 其次利用原型学习，其关键思想是基于与原型案例的相似性来解决新输入的问题，以检测共同点因果结构。 由此推断出$d-1$外生变量与目标变量之间的格兰杰因果关系和强度。

(Löwe 等人, 2022)提出了一种称为摊销因果发现(ACD)的框架，旨在训练单个模型来推断具有不同底层的样本之间的因果关系因果图但共享动态。 它是一个编码器-解码器框架，其中编码器函数被定义为推断输入样本的格兰杰因果关系，而解码器函数则根据推断的因果关系学习预测下一个时间步。 在实现中，图神经网络应用于摊销编码器，ACD使用变分推理对函数进行建模，该模型基于广泛使用的神经关系推理（NRI）模型（Kipf等人，2018）. 此外，为了导出推断边缘的因果解释，ACD 中提供了将零边缘函数与格兰杰因果关系联系起来的证明。 因此，无需重新拟合模型，就可以推断出先前未见过的样本的因果关系。

3.5. 其他的

上述四类方法一直是因果发现研究中许多努力的主题。 为了完整起见，我们在本小节中提出了与上述方法不同的五种方法，包括基于信息论统计的因果关系、基于微分方程的因果模型、非线性状态空间方法、基于逻辑的方法、和混合方法。

3.5.1. 基于信息论统计的因果关系

MTS 中​​的因果关系可以基于信息论统计来测量。 作为一种无模型度量，它广泛应用于基于约束的方法 (3.1.1) 和 Granger 因果模型 (3.4.2)。 但其定义和特点尚未详细说明。 在这一部分中，我们将首先介绍传递熵（Schreiber，2000），它是因果关系信息论统计的原始概念，然后是它的变体。

传递熵（Schreiber，2000）是对两个进程之间的信息流或有效耦合的度量，无论实际的功能关系如何。 作为一种无模型测量，它可以与特定结构学习方法的变体相结合，而不是基于模型的标准，该标准存在模型可能被错误指定的问题。 具体来说，从$i$到$j$（有时间滞后）的传递熵可以表示为：

$$\mathrm{TE}(X_i^t\to X_j^{t+1})=h(X_j^{t+1}|X_j^t)-h(X_j^{t+1}|X_j^t,X_i^t),$$

其中$h(\cdot |\cdot )$表示条件熵。 这里，术语$h(X_j^{t+1}|X_j^t)$测量给定$X_j^t$信息的$X_j^{t+1}$的不确定性，$h(X_j^{t+1}|X_j^t,X_i^t)$测量$X_j^{t+1}$的不确定性> 给出有关 $X_j^t$ 和 $X_i^t$ 的信息。 因此，当$X_i$的当前动态加上$X_j$的未来动态的不确定性时，我们可以从因果关系的角度来理解转移熵$\mathrm{TE}(X_i^t\to X_j^{t+1})$。 $X_j$。 对于高斯变量，转移熵与格兰杰因果关系的等价性如(Barnett等人, 2009)所示。 此外，在（Runge等人，2012b）中，将传递熵重新表述为分解形式并嵌入到多元图模型框架中。 (Runge 等人, 2012a)中，多元时间序列的因果耦合强度是基于传递熵的变体来量化的。

尽管在多变量场景中得到了一些利用，但传递熵仍受到成对限制的影响。 据报道，它无法区分网络中的直接因果关系和间接因果关系（Sun等人，2015）。 作为对成对限制的补救措施，引入了因果熵（Sun and Bollt，2014），这是一种用于推断因果关系的无模型信息理论统计量。 具体来说，以节点集$C$为条件的从节点集$I$到节点集$J$的因果熵定义如下：

$$\mathrm{CE}(X_{𝐈}^t\to X_{𝐉}^{t+1}|X_{𝐂}^t)=h(X_{𝐉}^{t+1}|X_{𝐂}^t)-h(X_{𝐉}^{t+1}|X_{𝐂}^t,X_{𝐈}^t),$$

这里$I,J,C$都是节点$\{1,2,\mathrm{\dots },d\}$的子集。 作为条件互信息的一种，因果熵是转移熵的推广，用于测量多变量网络关系的成对关系。 与传递熵和格兰杰因果关系的等价关系类似，因果熵在应用于高斯变量时也概括了格兰杰因果关系和条件格兰杰因果关系。 然而，根据其定义，该概念假设隐藏动态遵循平稳一阶马尔可夫过程，因为因果熵仅对时滞等于 1 的因果关系进行建模。 最近，为了测量任何滞后或瞬时关系，（Assaad 等人，2022a）中提出了因果熵的扩展，称为贪婪因果熵。

3.5.2. 基于微分方程的因果模型

微分方程是许多领域常用的建模工具，如果可以在相关时间尺度上进行测量，则特别有用。 与上述因果模型相比，此类方法专门用于对可以用微分方程很好地表示的系统进行建模（Peters等人，2022）。 在这一部分中，我们将首先回顾离散时间和连续时间的微分方程和因果模型之间的关系。 将引入第一个基于差异的因果发现框架。 然后，我们将给出此类方法的最新进展。

有大量文献(Peters 等人, 2022; Bongers 等人, 2018) (Schölkopf, 2019; Mooij 等人, 2013; Rubenstein 等人, 2018) 讨论微分方程和结构因果模型之间的关系。 对于离散时间，（Voortman等人，2010）首次提出了基于差异的因果发现框架。 满足跨时间限制，其中跨时间的所有因果关系都是由于导数 $\dot{x}$ 导致其积分 $x$ 发生变化。 这一特性使得基于差异的因果模型成为动态 SEM 的受限形式。 利用基于差异的因果学习器(DBCL)从数据中提取基于差异的因果模型，事实证明该模型能够识别反馈循环的存在或不存在。 对于连续时间，人们也做出了一些理论努力，通过常微分方程 (ODE) 推导出动态系统的因果解释（Mooij 等人, 2013; Rubenstein 等人, 2018; Blom 等人, 2019; Pfister等人，2019）和随机微分方程（SDE）（Hansen 和 Sokol，2014；Mogensen 等人，2018）。

最近，在多元时间序列不规则采样（在无穷小的时间间隔内）的动态因果系统下，（Bellot 等）提出了一种称为神经图形模型(NGM)的算法人，2022）。 在许多应用中，感兴趣的底层因果系统可以表示为动态结构模型，如下所示：

$$d𝐱(t)=𝐟(𝐱(t))dt+d𝐰(t),𝐱(0)={𝐱}_0,t\in [0,T],$$

其中$𝐰(t)$是$d$维标准布朗运动，${𝐱}_0$是独立于$𝐰(t)$的高斯随机变量，函数$𝐟$ 描述了因果图$G$。 NGM 是一种基于惩罚神经常微分方程 (neural-ODE) 的学习算法。 因果图的恢复可以转化为以下形式的惩罚优化问题：

$${\min }_{{𝐟}_{\theta }}\frac{1}{n}\sum _{i=1}^n{\Vert 𝐱(t_i)-\widehat{𝐱}(t_i)\Vert }_2^2,\mathrm{subject}\mathrm{to}{\rho }_{n,T}({𝐟}_{\theta })\mathrm{and}\widehat{𝐱}(t)={𝐟}_{\theta }(\widehat{𝐱}(t_i))dt,$$

其中系统的观察是在不规则的时间点。

3.5.3. 非线性状态空间方法

在这一部分中，我们将首先介绍非线性状态空间方法的基础知识，包括Takens定理和收敛交叉映射算法。 然后将给出原始算法的变体和最新进展，以解决诸如对噪声高敏感度、大样本需求、结果不一致和错误识别等挑战。

Takens提出的状态空间重构理论（Takens，1981）为分析非线性系统的动态特性提供了理论基础。 基于这一理论，另一种确定因果关系的方法，称为收敛交叉映射(CCM)，首先在（Sugihara等人，2012）中提出。 该方法专为耦合时间序列而开发，通过状态空间重建利用 Takens 定理。 具体来说，在给定两个时间序列 $x_1^t$ 和 $x_2^t$ 的情况下，首先分别使用 $x_1^t$ 和 $x_2^t$ 重建吸引流形 ${ℳ}_{x_1},{ℳ}_{x_2}$ 。 其次，可以通过测量${ℳ}_{x_1}$和${ℳ}_{x_2}$之间的对应关系来检测因果关系，具体来说，通过测试在一个流形上定义的每个局部邻域是否在另一个流形中保留。 图7给出了CCM的说明。 该方法已成功应用于许多非线性系统动态耦合的领域（Hirata等人，2016；Ye等人，2015）。

然而，原始CCM方法存在对观测噪声敏感度高、需要较多观测数据、不同优化算法结果不一致等问题。 为了克服这些挑战，基于时滞分析（叶等人，2015）、深度高斯过程（冯等人，2019）储层计算(Huang 等人, 2020b) 和神经 ODE (Brouwer 等人, 2021) 被提出。 此外，大多数基于 CCM 的方法最初是为双变量分析而开发的。 尽管可以多次使用相同的过程来确定多元时间序列之间的因果网络，但在高维条件下无法保证性能（黄等人，2020a）。 将间接因果关系误认为是直接因果关系是多变量环境中的关键挑战之一。 最近，提出了将 CCM 与偏相关相结合的部分交叉映射（PCM）（Leng 等人，2020）来消除间接因果影响。

图 7. 交叉收敛映射过程的图示。

3.5.4. 基于逻辑的方法

另一种方法用于时间序列数据中的因果推理和因果发现，它基于逻辑公式。 首先将介绍此类方法的原始算法，并结合其符号学和潜在因果关系的定义。 然后我们将给出它的变体和最新进展。

在基于逻辑的方法中，时态数据可以被视为对系统已占据的状态序列的观察，并在模型检查中被称为踪迹。 这一系列研究源于（Kleinberg 和 Mishra，2009） 中的工作，其中因果关系是用时序逻辑公式描述的。 具体来说，它首先利用逻辑，概率计算树逻辑(PCTL)，基于时间优先级和条件概率的提升来定义初步（潜在）因果关系。 给定原始作品中的符号，如果以下条件全部成立，则定义初步原因：(1) $F_{>0}^{\le \mathrm{\infty }}c$、(2) $c{\rightsquigarrow }_{\ge p}^{\ge 1,\le \mathrm{\infty }}e$ 和 (3) ，意味着$c$和$e$之间可能存在任意数量的转换，并且一组路径概率之和至少为$p$。为了将潜在的表面（潜在）原因分为真实原因和虚假原因，通过计算每个表面原因的平均概率差异引入了$ϵ$-不重要原因的概念与该效果的所有其他表面原因相关的效果：

$${ϵ}_{avg}(c,e)=\frac{{\sum }_{x\in X\mathrm{\backslash }c}{ϵ}_x(c,e)}{|X|},$$

其中${ϵ}_x(c,e)=P(e|c\wedge x)-P(e|\neg c\wedge x)$。 如果，表面原因$c$是$e$的$ϵ$-不重要的原因。 $ϵ$ 的值是根据经验原假设检验选择的，假设：(1) 数据包含两个类别：显着类别和不显着类别，(2) 显着类别相对于不显着类别而言较小。 同时实施错误发现率控制。 该方法也适用于领域（Kleinberg，2013）。

为了将方法扩展到离散和连续分量同时存在的条件，（Kleinberg，2011）中引入了PCTLc来表达涉及离散和连续变量的时间和概率属性，以及关系中的重要性。连续情况是通过效果的条件期望而不是条件概率来验证的。 此外，还提出了这种基于逻辑的方法的变体（Huang and Kleinberg，2015），以提高因果发现的准确性，并通过显示可以降低计算复杂度来更快地计算因果意义。几个条件。 遵循这种时序逻辑形式，最近的一项工作（da Costa and Dasgupta，2021）结合了决策树的思想，重新考虑因果发现问题，从实时时间中提取时序因果序列关系系列。

3.5.5. 混合方法：结合基于分数和基于约束的方法

混合方法的提出是为了结合基于约束的 (3.1) 和基于分数的 (3.2) 方法的优点。 我们涵盖了混合方法的两部分，包括基于最大-最小爬山启发式的方法，以及结合条件独立性测试来改进局部搜索的方法。

一些研究人员开发了基于最大-最小爬山启发式的混合方法（Tsamardinos等人，2006；Li和Ngom，2013；Li等人，2016）。 作为混合局部学习方法，最大-最小方法融合了基于约束的技术来限制潜在结构的空间和搜索和评分贝叶斯方法来搜索最佳结构的概念。 它们最初被用于静态数据的 BN 结构学习（Tsamardinos 等人，2006）。 （李等人，2016）中介绍的最大-最小爬山贝叶斯网络（MMHO-DBN），基于最大-最小爬山启发式的扩展来学习 DBN 的结构并用于真实基因表达时间序列数据的建模。

还有一些混合方法，将条件独立测试和局部搜索相结合来提高标准分数（Ogarrio等人，2016；Malinsky和Spirtes，2018）。 Greedy FCI (GFCI) (Ogarrio 等人, 2016) 是一种结合了 GES 和 FCI 特征的混合评分。 SVAR-GFCI （Malinsky 和 ​​Spirtes，2018） 将此方法扩展到时间序列的因果结构学习。 在（Sanchez-Romero等人，2019）中，提出了一种称为快速邻接偏度（FASK）的PC稳定算法的变体，以及一种混合两步算法来提取因果关系对于时间序列数据。

4. 从事件序列中发现因果关系

多元时间序列中的一个重要假设是时间戳是离散的并且时间间隔是固定的。 然而，在现实世界中，绝大多数事件不会以固定的时间间隔发生。 因此，我们需要想出一些方法来处理这些不规则的、异步的数据。 我们可以将事件序列构造为$\{(t_1,e_1),(t_2,e_2),\mathrm{\dots }\}$，其中第一维表示对应事件发生的时间，第二维表示对应事件类型。 在本节中，我们将重点关注推断事件序列中的因果关系。 首先，介绍多元点过程，这是事件序列中因果发现的基础知识。 然后，我们回顾基于格兰杰因果模型的方法，这些方法已经很成熟。 最后，给出了其他方法，包括基于约束和基于分数的方法。

4.1. 多元点过程

事件序列记录一种特定类型事件（或简单起见，称为“事件类型”）的发生。 同时，我们可以通过点过程来表征事件序列。 为了发现不同类型事件之间的关系，我们考虑其高维情况，即通过多元点过程（MPP）对事件序列进行建模。 因此，我们的问题可以定义为输入一组点过程，其中每个点过程代表一个事件序列，并输出由不同过程建立的因果图$G$。 在因果图中$G$，每个节点代表一个点过程，每个有向边捕获从一个点过程到另一个点过程的有向交互。 在这一部分中，我们将详细介绍 MPP，包括它们的强度函数和对数似然函数。

MPP 的强度函数。 时间点过程是由连续时间内发生的二元事件的时间序列组成的随机过程（Daley等人，2003）。 MPP 是高维点过程，这意味着它们可以涉及多种类型的事件。 $ℰ=\{1,\mathrm{\dots },E\}$ 是事件类型的集合。 这些事件$\{t_1,t_2,\mathrm{\dots },t_n|t_i\in [0,T]\}$的发生时间分布不均匀。 具有$E$类型事件的多元点过程可以用$E$计数过程${\{N_e\}}_{e=1}^E$来表示，其中$N_e=\{N_e(t)|t\in [0,T]\}$。 点过程的核心是其条件强度函数，其中捕获了过程的模式。 type-$u$ 强度函数可以定义为给定历史的 type-$e$ 事件发生的预期瞬时速率：

$${\lambda }_e(t)=\frac{𝔼[d{N}_e(t)|{\mathscr{H}}_t]}{dt}$$

这里${\mathscr{H}}_t=\{(t_i,e_i)|t_i<t,e_i\in ℰ\}$表示在时间$t$之前发生的所有类型的事件。

MPP 的对数似然函数。 接下来，我们展示联合分布的强度函数和概率密度函数（PDF）之间的关系：$f((t_1,e_1),\mathrm{\dots },(t_n,e_n)|(t_0,e_0))$。 使用链式法则，有$f((t_1,e_1),\mathrm{\dots },(t_n,e_n)|(t_0,e_0))={\prod }_{j=1}^{n}f((t_j,e_j)|{\mathscr{H}}_{t_j})$。 然后，我们可以设置似然函数来估计联合分布：

(1) $$\widetilde{L_0}\triangleq \sum _{j=1}^{n}lnf(t_j|e_j,{\mathscr{H}}_{t_j})+\sum _{j=1}^{n}lnf(e_j|{\mathscr{H}}_{t_j})$$

由于目标是推断不同事件之间的因果关系，因此这里我们专注于第一项，省略第二项：$L_0\triangleq {\sum }_{j=1}^{n}lnf(t_j|e_j,{\mathscr{H}}_{t_j})$。 强度函数反映了在给定${\mathscr{H}}_{n-1}\triangleq {\mathscr{H}}_{t_{n-1}}$信息的情况下对$[t,t+dt]$中发生的事件的期望。 与生存分析中死亡率的计算类似，有：

${\lambda }_e(t|{\mathscr{H}}_{n-1})$ $={\displaystyle \frac{𝔼[N_e(t+dt)-N_e(t)|{\mathscr{H}}_{n-1}]}{dt}}$

$={\displaystyle \frac{{\int }_{t_{n-1}}^{t+dt}f(l|e,{\mathscr{H}}_{n-1})𝑑l-{\int }_{t_{n-1}}^{t}f(l|e,{\mathscr{H}}_{n-1})𝑑l}{1-{\int }_{t_{n-1}}^{t}f(l|e,{\mathscr{H}}_{n-1})𝑑l}}$

$=-{\displaystyle \frac{d}{dt}}ln(1-{\displaystyle {\int }_{t_{n-1}}^t}f(l|e,{\mathscr{H}}_{n-1})𝑑l)$

对上面的方程进行积分，并将$f$的表达式代入1，我们有，

(2) $$L_0=\sum _{j=1}^n\{ln{\lambda }_{e_j}(t_j|{\mathscr{H}}_{j-1})-{\int }_{t_{j-1}}^t{\lambda }_{e_j}(l|{\mathscr{H}}_{j-1})𝑑l\}$$

我们在上面简要介绍了多元点过程并构造了似然函数来表征 MPP。 接下来，我们的目标是使用基于 Granger 的方法以及基于约束和基于分数的方法来发现 MPP 中的因果关系。

4.2. 基于格兰杰因果关系的方法

在本小节中，我们考虑推断事件序列中格兰杰因果关系的任务。 与 MTS 中​​的情况类似，如果 对于预测 $e_i(t)$ 有用，我们称 $e_j$ 类型事件 Granger 导致 $e_i$。 具体方法可根据以下模型规格进行分类：即、GLM 点过程、霍克斯过程、沃尔德过程和神经点过程。

4.2.1. GLM 点过程的方法

我们首先介绍通过点过程 (Truccolo et al, 2005) 的 Generalized Linear Model (GLM) 建模的事件序列因果发现方法。 GLM 假设强度函数的对数具有线性格式，即$ln{\lambda }_e(t|{\mathscr{H}}_{n-1})={\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2$。 具体来说，在我们的任务中，强度函数如下：

(3) $$ln{\lambda }_i(t|{\gamma }_i,H_i(t))={\gamma }_{i,0}+\sum _{j=1}^J\sum _{m=1}^{M_i}{\gamma }_{i,j,m}{R}_{j,m}(t)$$

这里${\gamma }_{i,0}$可以解释为事件$e_i$的背景强度，${\gamma }_{i,j,m}$是由type-触发的type-$e_i$事件的强度$e_j$ 事件，$R_{j,m}(t)$ 是 $[t-mW,t-(m-1)W]$ 中发生的 $e_j$ 类型事件的发生次数($W$ 是一个小数字，指的是时间范围的长度）。 通过观察${\sum }_{m=1}^{M_i}{\gamma }_{i,j,m}$的符号，我们可以区分$e_j$类型事件对$e_i$类型事件是否具有兴奋或抑制作用。

为了推断类型 $e_j$ 和类型 $e_i$ 事件之间的格兰杰因果关系，我们将 3 代入似然函数 2。 接下来，我们遵循一个简单的想法，即我们可以排除某种类型的事件，然后通过将其强度与原始情况进行比较来推断格兰杰因果关系。 具体来说，我们获得了在有和没有 type-$e_j$ 影响的情况下 type-$e_i$ 出现的可能性：$L_i({\gamma }_i)$、$L_i({\gamma }_i^j)$。 然后，考虑${\sum }_{m=1}^{M_i}{\gamma }_{i,j,m}$是影响类型的指示符，从类型-$e_j$到类型-$e_i$事件的格兰杰因果关系可以建议为 （Kim等人，2011）：

${\varphi }_{ij}$ $=-sign({\displaystyle \sum _{m=1}^{M_i}}{\gamma }_{i,j,m}){\mathrm{\Gamma }}_{ij}$

$=-sign({\displaystyle \sum _{m=1}^{M_i}}{\gamma }_{i,j,m})log{\displaystyle \frac{L_i({\gamma }_i^j)}{L_i({\gamma }_i)}}$

显然，存在$L_i({\gamma }_i)\ge L_i({\gamma }_i^j)$，因此存在${\mathrm{\Gamma }}_{ij}=log\frac{L_i({\gamma }_i^j)}{L_i({\gamma }_i)}\le 0$。 仅当满足“”时，类型$e_j$事件才会成为类型$e_i$事件的格兰杰原因。 接下来，Kim 等人 (Kim 等人, 2011) 通过进行 $H_0$ 假设：${\theta }_0={\gamma }_i^j$ 和 $H_1$假设：${\theta }_1={\gamma }_i$。 通过 FDR 显著性检验，最终的因果关系可以用 $\widetilde{{\varphi }_{ij}}$ 来估计：（1）当 $\widetilde{{\varphi }_{ij}}>0$ 时，类型-$e_j$ 事件是类型-$e_i$ 事件的兴奋原因；（2）当 时，原因是抑制性的；（3）当 $\widetilde{{\varphi }_{ij}}=0$ 时，类型-$e_j$ 和类型-$e_i$ 事件之间不存在因果关系。

4.2.2. 霍克斯过程的方法

在这一部分中，我们回顾霍克斯过程的方法。 作为一种特殊类型的点过程，首先给出霍克斯过程的基础知识。 然后我们详细介绍基于 MLE 推断因果关系的方法，包括（1）参数化策略和（2）正则化方法。 接下来，我们回顾其他估计方法，包括（1）图形事件模型，（2）广义矩方法，（3）事件序列分离，以及（4）最小描述长度。 我们注意到，由于格兰杰因果关系和霍克斯过程之间的自然匹配，这一类别中存在大量文献。

霍克斯过程是一种点过程，具有固定形式的强度函数：

(4) $${\lambda }_{e_i}(t)={\mu }_{e_i}+\sum _{e_j=1}^E{\int }_0^t{\varphi }_{e_i{e}_j}(s)𝑑N_{e_j}(t-s)$$

这里，${\mu }_{e_i}$被称为基线强度，它只能受到外源事件的影响，因此，随着时间的推移是一个常数。 影响函数${\varphi }_{e_i{e}_j}(s)$衡量由历史类型$e_j$事件触发的未来类型$e_i$事件的兴奋度的衰减。 也就是说，它捕获从$e_j$到$e_i$的内源强度。 考虑到$\varphi$的定义与格兰杰因果关系的相似性，我们可以通过分析$\varphi$直接推断出格兰杰因果关系：

命题1。

(Eichler 等人 (Eichler 等人, 2017), 2017)

$$e_j\text{does not Granger-cause}{e}_i\iff {\varphi }_{e_i{e}_j}(s)=0,\forall s\in R$$

因此，我们的目标是为每个事件和所有$t\in R$建模${\varphi }_{e_i{e}_j}(t)$。 然而，由于事件序列的复杂性和异质性，这一任务可能极其难以完成。 周等人 (Zhou 等人, 2013) 将 ${\varphi }_{e_i{e}_j}(s)$ 参数化为 $a_{e_i{e}_j}g(s)$。 通过这种方式，$\varphi$被分成事件交互部分和时间衰减部分。

MLE 方法。 可以执行最大似然估计（MLE）来估计4中的参数。 我们将$\lambda$的表达式代入2，得到相应的似然函数：$L(A,\mu )$。 这里，$A$由$(a_{e_i{e}_j})$组成，$\mu$由${\mu }_{e_i}$构建。 接下来，考虑到在现实场景中，大多数事件只能影响一小部分其他事件，并且影响网络中的社区结构往往排名较低（周等人，2013） ，我们应该向 MLE 损失函数添加惩罚项。 具体来说，为了实现矩阵$A$的低秩和稀疏性，可以构造以下目标函数：

$$\underset{A\ge 0,\mu \ge 0}{min}-L(A,\mu )+{\lambda }_1{\Vert A\Vert }_{\ast }+{\lambda }_2{\Vert A\Vert }_1$$

这里，$\parallel \cdot {\parallel }_{\ast }$是核范数，其降低矩阵秩的性能已经被证明。 而$\parallel \cdot {\parallel }_1$就是L1范数。 它可以强制矩阵 A 获得更多稀疏性。 ${\lambda }_1$、${\lambda }_2$是控制这两个惩罚强度的参数。 我们将对象函数表示为$f(A,\mu )$。 除此之外，可以采用基于EM的算法来解决$A$和$\mu$的优化问题。 具体来说，周使用代理函数$Q(A,\mu ;A^{(m)},{\mu }^{(m)})$作为$f(A,\mu )$的紧上限。 通过迭代优化$Q(A,\mu ;A^{(m)},{\mu }^{(m)})$,$f(A,\mu )$被迫减少，从而成功优化。 我们在算法 LABEL:alg:Algorithm1 中总结了这一点。

（1）参数化策略：如果数据不符合其强假设，上面提到的参数化方法${\varphi }_{e_i{e}_j}(s)$ = $a_{e_i{e}_j}g(s)$可能会出现性能不佳的情况。 因此，为了在不同类型的事件序列中获得鲁棒性，Xu 等人 (Xu 等人, 2016) 提出了选择一系列基本函数并使用它们的线性组合的策略，${\sum }_{m=1}^M{a}_{e_i{e}_j}^m{\kappa }_m(s)$，对目标强度函数进行建模。

然而，在NPHC (Achab 等人, 2017)中，Achab 等人提出，基本函数策略在存在时将具有非凡的计算复杂度。存在太多类型的事件（即 $E$ 很大）。 由于我们的目标只是推断格兰杰因果关系，因此不需要完全参数化霍克斯过程，因此，我们只需要估计相应的积分${\int }_0^{+\mathrm{\infty }}{\varphi }_{e{e}^{\prime }}(s)𝑑s$。 Achab将积分表示为$g_{e{e}^{\prime }}$，而$(g_{e{e}^{\prime }})$构成矩阵G。那么，由Eichler的证明(Eichler等人, 2017)以及${\varphi }_{e{e}^{\prime }}(s)>0,\forall s>0$，显然是$(g_{e{e}^{\prime }})=0\iff e^{\prime }\text{does not Granger-cause}e$。

其他作品考虑了事件序列内的潜在拓扑关系。 在THP (Cai 等人, 2021)中，Cai 等人假设事件之间存在隐藏的无向图结构$G_N$。 相应的强度函数为${\lambda }_{e_i}(n,t)={\mu }_{e_i}+{\sum }_{e_j\in E}{(g_{e_i{e}_j}\ast s_{e_i,e_j,t})}_{G_N}(n)$。 这里，$g_{e_i{e}_j}$是图卷积核，可以捕获图邻居的影响。 而$s_{e_i,e_j,t}$是时间卷积核，代表过去的影响函数之和。这是基于隐藏的拓扑结构在此过程中不会改变的假设。

（2）正则化方法：在上述方法中，我们提出了$A$的核范数${\Vert A\Vert }_{\ast }$和L1范数${\Vert A\Vert }_1$作为正则化器。 在基本函数方法中，应用特殊的稀疏组套索正则化器（Simon等人，2013）来拟合其求和参数化。 具体来说，Xu 等人（Xu 等人，2016）同时进行了分组拉索惩罚和拉索惩罚，以便对所有$m$执行$a_{e_i{e}_j}^m=0$，即分组稀疏性，此外还对所有条目$a_{e_i{e}_j}^m$执行常规稀疏性。 然而，在 $L_0$ Hawkes (Idé 等人, 2021) 中，Ide 等人证明了具有 L1 正则化的基于 EM 的 MLE 算法无法在数学上提供稀疏解。 因此，它们的稀疏解只能显示为数值伪影。 随后，Ide 提出了 L0 正则化 EM-MLE 算法来解决这个问题。 这里，L0范数${\Vert A\Vert }_0$表示矩阵$A$中非零条目的数量。

与上述拓扑参数化策略类似，Xu 等人 (Xu 等人, 2016) 在构建我们的正则化器时考虑了事件类型之间的潜在拓扑关系。 具体来说，可以呈现成对相似性${\sum }_{e_i=1}^E{\sum }_{e_j\in C_e}{\Vert a_{e_i\cdot }-a_{e_j\cdot }\Vert }_F^2+{\Vert a_{\cdot e_i}-a_{\cdot e_j}\Vert }_F^2$，以强制相似的事件可以具有相似的强度函数。 然而，我们必须补充一点，这个正则化器需要预定义的集群结构，因此可以进行优化。

在发现格兰杰因果关系时，先前的领域知识可能非常有用。 由于事件序列的高维性和异质性，现有算法经常出现拟合不足和可解释性差的问题。 因此，很自然地考虑将人类的领域知识添加到因果推断模型中。 具体来说，建立了一个带有用户反馈的自下而上的可视化模型（Jin等人，2021）。 Jin等人在MLE-SGLP中用传统的MLE方法建立了他们的基础模型。 在训练过程中，用户可以根据网络中的领域知识确认或删除因果关系。 并且模型会根据用户的选择改变其优化目标。 例如，根据MLE-SGLP (Xu 等人, 2016)中的思想，Jin (Jin 等人, 2021)将其强度函数构造为${\varphi }_{e_i{e}_j}={\sum }_{m=1}^M{a}_{e_i{e}_j}^m{\kappa }_m(s)$，并将$a_{e_i{e}_j}$设置为$[a_{e_i{e}_j}^1,\mathrm{\dots },a_{e_i{e}_j}^n]$。 相应地，他们的目标函数可以是：$\underset{\mu ,\alpha }{argmin}-L+\alpha {\sum }_{e_i,e_j}{\Vert a_{e_i{e}_j}\Vert }_2$。 当用户选择确认或删除因果图中的边 $\widehat{G}$ 后，Jin 将对象函数更新如下：

(5) $$\begin{array}{c}\hfill \underset{\mu ,\alpha }{argmin}-L+{\alpha }_v\sum _{e_i,e_j}{\Vert a_{e_i{e}_j}(\widehat{G})\Vert }_2\\ \hfill \text{s.t.}{a}_{e_i{e}_j}=0\text{for}(e_j\to e_i)\notin \widehat{G}\end{array}$$

(6) $$a_{e_i{e}_j}(\widehat{G})=\{\begin{array}{cc}0;\text{if}(e_j\to e_i)\text{is confirmed}\hfill & \\ a_{e_i{e}_j};\text{otherwise}\hfill & \end{array}$$

这里5中的约束适合删除操作，6中的更新代表用户的确认。

其他估计方法

（1）图形事件模型：上述方法使用最大似然估计来对事件序列的霍克斯过程进行建模。 然而，这些尝试缺乏可解释性，需要对参数进行微调过程才能实现良好的性能。 因此，可以提出完全数据驱动、基于图形和依赖性捕获的图形事件模型（GEM）来推断事件序列中的格兰杰因果关系。

我们将在4.3中详细阐述GEM的属性。 这里，我们只关注它与格兰杰因果关系的关系。 假设有一个有向图$𝒢=(ℰ,𝒜)$，其中的边代表不同事件类型之间的依赖关系。 对于每个事件类型 $e$，我们假设其条件强度只能受其父类型影响，即它遵循 ${\lambda }_e(t|h_t)={\lambda }_e(t|{[h_t]}_{P_a(e)})$，其中 $P_a(e)\subseteq ℰ$ 是 $e$是图表$𝒢$中的父事件，${[h_t]}_{P_a(e)}$是集合$P_a(e)$中列出的类型的事件的历史记录。 根据1，有，

命题2。

（GEM中的格兰杰因果关系，Yu 等人，2020 （Yu 等人，2020）)

对于两种事件类型$e_i$ $𝒢=(ℰ,𝒜)$ 中的 > 和 $e_j$,$e_j$ 不会导致 $e_i\iff {\varphi }_{e_i{e}_j}(t)=0,\forall t>0\iff e_j\notin Pa(e_i)$ 的格兰杰原因

因此，人们可以应用传统的基于分数的结构学习方法来发现格兰杰因果关系。 例如，可以呈现BIC分数以用于学习优化图${𝒢}^{\ast }$。 优化方法是一致的。 同时，Yu 进行了前向后向搜索（FBS）来独立学习某种事件类型的父类型（Yu 等人，2020）。 对于 GEM 家族（Gunawardana 和 Meek，2016），使用 BIC 分数的前向-后向搜索被证明是合理且完整的。

（2）广义矩量法：在NPHC中，优化对象是矩阵$G=(g_{e{e}^{\prime }})=({\int }_0^{+\mathrm{\infty }}{\varphi }_{e{e}^{\prime }}(s)𝑑s)$。 因此，广义矩法（GMM）可以用来解决这个问题（Hall，2004）。 Achab 等人提出了一种基于 GMM 的 NPHC 算法，对矩阵 $G$ 的一阶、二阶和三阶累积量进行建模（Achab 等人, 2017）。 之后，从$G$就可以直接得到格兰杰因果关系。这种矩估计方法被证明对于某些观测噪声是一致且稳健的（Trouleau等人，2021）。 然而，这种方法在特定数据集（例如长尾数据集）中可能会收到较差的结果。 这主要是由于 GMM 的普遍问题：它们只能捕获统计分布矩内的信息。

（3）事件序列分离：另一个有趣的想法是将事件序列分离为多个子序列，并在每个子序列中相应地应用霍克斯过程模型。 在GC-nsHP (Chen 等人, 2022a)中，Chen 等人将事件序列${\mathscr{H}}_n=[(t_1,e_1),\mathrm{\dots },(t_n,e_n)]$分为$K$不同模式，其中$K$需要根据其应用场景进行预定义。 相同模式中的“事件”应该构建 ${\mathscr{H}}_n$ 的固定子流程。 然后，专门针对$K$模式建立了$K$不同的霍克斯过程，并且只能在每个模式内部学习格兰杰因果关系。 在每次迭代中，交替执行基于维特比路径的模式重新分配算法和基于 EM-MLE 的参数更新算法。 在参数更新部分，考虑到$X_{t-1}$和$X_t$更有可能处于同一模式，Chen添加了惩罚项以帮助将相邻序列放入同一模式。

（4）最小描述长度：遵循最小描述长度（MDL）原则（Rissanen，1998；Grünwald and Roos，2019），Jalaldoust等人进行了权衡拟合优度和模型复杂度之间的关系（Jalaldoust 等人，2022）。 具体来说，他们将参数空间$\mathrm{\Theta }$划分为$\{{\mathrm{\Theta }}_{\gamma }:\gamma \in \mathrm{\Gamma }\}$，定义了幸运函数$v:\mathrm{\Theta }\to ℝ$，并为每个模型设置归一化最大似然分布$\gamma \in \mathrm{\Gamma }$ 为：

(7) $$p_{v|\gamma }^{NML}(x)=\frac{ma{x}_{\theta \in {\mathrm{\Theta }}_{\gamma }}p(x|\theta )v(\theta )}{{\int }_{x\in 𝒳}ma{x}_{\theta \in {\mathrm{\Theta }}_{\gamma }}p(s|\theta )v(\theta )𝑑s}$$

积分的对数可以看作模型复杂度：

(8) $$COMP(M_{\gamma };v)=log{\int }_{x\in 𝒳}ma{x}_{\theta \in {\mathrm{\Theta }}_{\gamma }}p(s|\theta )v(\theta )𝑑s$$

Jalaldoust 使用以下方法选择了优化模型 ${\widehat{\gamma }}^{MDL}\in \mathrm{\Gamma }$

(9) $${\widehat{\gamma }}^{MDL}=\underset{\gamma \in \mathrm{\Gamma }}{\arg \min }{L}_v(\gamma ;x)=\underset{\gamma \in \mathrm{\Gamma }}{\arg \min }(-log\pi (\gamma )-r_v({\widehat{\theta }}_{v|\gamma }(x);x)+COMP(M_{\gamma };v))$$

其中 $\pi$ 是均匀分布，$r_v({\widehat{\theta }}_{v|\gamma }(x);x)$ 是与 $p$、$v$ 和 ${\mathrm{\Theta }}_{\gamma }$

此外，考虑在所有二元 $p\times p$ 矩阵集合内从 $\gamma \in \mathrm{\Gamma }$ 到因果图的 $p\times p$ 相邻矩阵的一对一映射。 通过优化9，人们可以从预定义的模型族中选择最合适的模型，从而推断事件类型之间的格兰杰因果关系。

4.2.3. Wold 工艺方法

虽然大多数涉及从事件序列中发现格兰杰因果关系的现有算法都基于霍克斯过程，但我们也可以在另一种类型的过程 - Wold 过程上对这些关系进行建模，后者本质上复杂性较低。 假设我们将 ${\delta }_i=t_i-t_{i-1}$ 表示为从第 $(i-1)$ 个事件发生开始到第 $i$ 个事件的等待时间。 Wold 进程建立在一个简单的假设之上，即当前等待时间 ${\delta }_i$ 仅与最近的过去等待时间 ${\delta }_{i-1}$ 相关。 也就是说，集合$\{{\delta }_i,i\in \mathrm{N}\}$构成了一条马尔可夫链。 Wold 过程中固有的马尔可夫特性使其适合对某些复杂系统的动力学进行建模。 此外，Figueiredo 等人 (de Figueiredo 等人, 2018)) 在某些数据集上测量了 ${\delta }_i$ 和 ${\delta }_{i-1}$ 之间的相关性。 结果表明，在大多数数据集中，皮尔逊相关性中值均高于 0.7，这是 Wold 模型充分性的标志。 据此，根据 Alve 等人 (da Silva Alves 等人, 2016) 和 Figueiredo 等人 (de Figueiredo 等人, 2018) 的思想，强度函数可以为执行为

$${\lambda }_{e_i}(t)={\mu }_{e_i}+\sum _{e_j\in E}\frac{{\alpha }_{e_i{e}_j}}{{\beta }_{e_j}+{\mathrm{\Delta }}_{e_i{e}_j}(t)}$$

基于BuSca模型。 这里，${\mu }_{e_i}$ 是霍克斯过程中的基本强度。 ${\mathrm{\Delta }}_{e_i{e}_j}(t)$ 表示在时间 t 上最后一次发生 $e_i$ 类型事件与 $e_j$ 类型事件之间的时间间隔。也就是说，如果我们定义在时间 $t$ 之前最近的 $e_i$ 类型事件发生在时间 $t_{e_i}$，在时间 $t_{e_i}$ 之前最近的 $e_j$ 类型事件发生在时间 $t_{e_j}$，则有 ${\mathrm{\Delta }}_{e_i{e}_j}(t)=t_{e_i}-t_{e_j}$。 因此，如果 ${\mathrm{\Delta }}_{e_i{e}_j}(t)$ 减少，强度函数中的交叉类型条目 ${\sum }_{e_j\in E}\frac{{\alpha }_{e_i{e}_j}}{{\beta }_{e_j}+{\mathrm{\Delta }}_{e_i{e}_j}(t)}$ 将会更大。 这完全符合这样一个事实：如果 type-$e_i$ 事件总是在 type-$e_j$ 事件发生之前发生，我们就会看到 $e_j$ 具有一定效果的可能性更大在$e_i$上。 ${\alpha }_{e_i{e}_j}$是满足${\sum }_{e_j\in E}{\alpha }_{e_i{e}_j}=1$的规范化条目，而${\beta }_{e_j}$是基本速率，使得当两种类型之间的时间间隔${\mathrm{\Delta }}_{e_i{e}_j}(t)$无穷小时在时间$t$，交叉类型条目将收敛到$\frac{{\alpha }_{e_i{e}_j}}{{\beta }_{e_j}}$。

格兰杰因果关系可以通过检查 ${\alpha }_{e_i{e}_j}$ 通过这个基于 Wold 的模型来了解。 具体地，如果${\alpha }_{e_i{e}_j}\ne 0$，则认为$e_j$格兰杰原因$e_i$。 由于学习过程的方法可能没有足够稀疏的解决方案，菲格雷多测试了这些可能的格兰杰因果关系的统计显着性，并丢弃那些显着性较低的。 此外，基于Wold的模型可以通过MCMC、期望最大化（EM）（de Figueiredo等人，2018）和变分推理（Etesami等人，2021）来学习接近。 任务是推断强度函数中的参数$\{{\alpha }_{e_i{e}_j},{\beta }_{e_j},{\mu }_{e_i}|\forall e_i,e_j\}$，这可以揭示事件序列内的所有属性。 这里我们不详细阐述这些学习方法的细节。

4.2.4. 神经点处理方法

随着神经网络的快速发展，神经点过程（NPP）逐渐被用来建模事件序列和推断因果关系。 这些NPP算法的核心思想是利用神经网络来推断强度函数${\lambda }_e(t)$。 具体来说，它们将事件序列编码到隐藏状态，在此期间它们捕获序列的特征。 然后，他们使用解码器来推断未来的强度函数。 核电站有两种主要类型。 一种是基于自回归（AR）模型；它的隐藏状态 $h_i$ 仅在事件发生时更新。 另一种遵循隐藏状态$h(t)$随时间连续变化的假设。 连续时间模型的优点是自然，并且由于其连续的特性，更适合估计任何时间$t$的属性。 尽管如此，这种灵活性是有代价的。 与基于 AR 的离散模型相比，连续隐藏模型的训练速度可能会较慢。 这是因为演化、可能性和采样过程可能需要数值近似。 在这一部分中，我们将首先介绍 NPP 的基础知识，然后介绍如何学习 NPP 中的格兰杰因果关系。

核电站基础知识。 使用基于AR的NPP模型推断格兰杰因果关系（张等人，2020）的一般流程如下。 首先，我们将每个事件嵌入到一个向量$v_i=[\theta (t_i-t_{i_1});V^T{z}_i$]中，其中$\theta (\cdot )$是预定义函数，V是事件类型的嵌入函数，$[\cdot ;\cdot ]$可以是连接和 $z_i$ 可能是事件类型 $u_i$ 的 one-hot 编码。 然后，我们利用序列编码器（例如 LSTM 或 GRU）将 $\{v_j;j\le i\}$ 编码为 $h_i=Enc(h_{i-1},v_i)$。 此外，还存在一种不同的编码方法，其中使用例如自注意力（self-attention）等对每个 i 独立地进行编码。 这可以更好地捕获事件之间的远程依赖关系，但计算复杂度也很高。

接下来，我们的目标是将隐藏状态 $h_i$ 解码为强度函数 ${\lambda }_e(t)$。 为此，我们需要对 $\lambda$ 做出一些假设。 例如，我们可以预测强度函数可以分为一些与交互相关和与时间相关的函数的总和${\lambda }_e(t)={\sum }_{m=1}^M{a}_{em}{\kappa }_m(t)$，类似于Xu的做法（Xu等人，2016）. 接下来，我们只能推断 $a_{em}$ 条目，因为 ${\kappa }_m(s)$ 可以从一个大的函数族中选择，其中的函数可以表示各种时变模式。 因此，

$$\alpha :{ℝ}^{rank(h_i)}\to {ℝ}_{+}^{K\ast S}$$

是该模型对应的解码器。 其中，$k$和$S$分别为事件类型数量和相应的基本功能。 然而，上述方法无法拟合连续隐藏状态模型。 在这种情况下，由于$h(t)$是连续的，因此比$h_i$携带了更多的时变信息，我们可以简单地将强度定义为：

(10) $${\lambda }_e(t)=g_e(h(t))$$

这里，$g_e:{ℝ}^{rank(h_i)}\to R_{>0}$ 是一个非线性函数（例如软加函数），它将 $h(t)$ 映射到事件类型 $e$ 在时间 $t$ 时的相应强度函数。

至于训练过程，现在大多数NPP都像大多数传统的PP方法一样使用最大似然估计（MLE）方法。 他们以最大似然估计的负对数似然为目标函数，并利用神经网络对其进行优化。 此外，还有其他方法可用于学习 MLE。 例如，如果我们将目标设置为${𝔼}_{X\sim p(X)}[f(X)]$，我们可以通过变分推理或强化学习使用$f(X)$对点过程$P(X)$进行建模。

从 NPP 推断格兰杰因果关系。 在推断格兰杰因果关系时，需要使用归因方法（张等人，2020）。 这是因为，在神经方法中，大多数算法不遵循霍克斯过程中的参数化。 相反，他们的目标是直接对过程的强度函数进行建模，以放松霍克斯过程的严格性，从而获得更高的准确性。 由于这些强度函数捕获了事件序列的所有特征，因此我们应该充分利用它们。 为此，张等人首先将$x_p=[t_1,e_1,\mathrm{\dots },t_p,e_p,t_{p+1}]$、$\underset{¯}{x_p}=[t_1,0,\mathrm{\dots },t_p,0,t_{p+1}]$表示为基线输入，将$f_k(x_p)={\int }_{t_p}^{t_{p+1}}{\lambda }_e(s)𝑑s$表示为影响函数(Zhang 等人, 2020 ）。 对于每个事件类型 k，我们有：

$$f_k(x_p)-f_k(\underset{¯}{x_p})=\sum _{q=1}^p{A}_q(f_k,x_p,\underset{¯}{x_p})$$

其中 $A_q(f_k,x_p,\underset{¯}{x_p})$ 是归因（例如 积分梯度）用于事件类型 $z_q$。 因此，$A_q(f_k,x_p,\underset{¯}{x_p})$可以被视为$z_j$类型事件对给定历史$x_p$的k类型事件预测的贡献。 接下来，张对$A_q(f_k,x_p,\underset{¯}{x_p})$进行归一化为

$$Y_{e_i,e_j}=\frac{{\sum }_{s=1}^s{\sum }_{p=1}^{n_s}{\sum }_{q=1}^i𝕀(k_q^s=e_j)A_q(f_{e_i},x_p^s,\underset{¯}{x_p^s})}{{\sum }_{s=1}^s{\sum }_{p=1}^{n_s}𝕀(k_q^s=e_j)}$$

接下来，可以从$Y_{e_i{e}_j}$推断出$u_i$和$u_j$类型事件之间的格兰杰因果关系。 该方法不仅可以测量抑制因果关系，还可以测量因果关系的大小。

有趣的是，其他一些神经算法只是像霍克斯过程中那样对强度函数进行建模。 他们将$\mu$和$\alpha$设置为矩阵，并将$H$和$A$直接放入神经网络中。 由于输入结构更容易，我们可以添加其他假设，例如事件之间的拓扑结构，并让神经网络（在本例中为 GCN）迭代优化 $H$ 和 $A$ 。 然后，我们可以直接从矩阵$A$推断格兰杰因果关系。

4.3. 其他推理方法

在本节中，我们不会直接对点过程中的强度函数进行建模。 相反，我们专注于发现不同流程（即不同类型的事件）之间的关系。 为此，我们可以利用前面提到的图形事件模型，并放松每个节点遵循霍克斯流程的假设。 从历史上看，Didelez 等人和 Meek 等人首先引入了图事件模型来捕获事件之间的依赖关系。 基于常见的图方法，他们假设事件类型的强度函数仅与其父类型相关。 GEM 随着时间的推移捕获各种类型事件之间的依赖关系，提供了一个通用框架来对图方法中的依赖关系进行建模。 因此，与平稳和离散时间情况类似，可以利用基于约束和基于分数的方法。

4.3.1. 基于约束的方法

就像不同随机变量之间的独立性概念一样，我们可以为点过程定义过程独立性：

定义1。

(Didelez, 2008 (Didelez, 2008); Bhattacharjya 等人, 2022 (Bhattacharjya 等人, 2022))

对于进程 $X,Y,Z$ ， 英石。 如果 $X$ 中的所有事件都具有条件​​强度，则 $Y\cap Z=0$、$X$ 是一个独立于给定 $Z$ 的 $Y$ 的过程如果 $Z$ 中事件的历史信息已知，则进程 $Y$ 中的事件不会提供任何进一步的信息。

Meek 等人 (Meek, 2014) 和 Bhattacharjya 等人 (Bhattacharjya 等人, 2022) 引入了 ${\delta }^{\ast }$-分离的概念，它基于on $d$-separation 但释放了不具有自循环的限制，并使每个自循环独立于其自己的历史。 然后，他们提出了类似于忠实度假设的 ${\delta }^{\ast }$ 分离的因果依赖假设。 基于因果依赖假设，提出了几种基于约束的方法，例如PC和最大最小父母算法，来学习不同类型之间的因果关系。 有多种过程独立性测试仪可供选择。 例如，我们有 NI 测试仪：

(11) $$\text{NI score}=\frac{1}{2}\frac{{\sum }_z{({\lambda }_{x|y,z}-{\lambda }_{x|\widehat{y},z})}^2}{{\sum }_z{\lambda }_{x|y,z}+{\lambda }_{x|\widehat{y},z}){}^2}$$

其中 $y$ 和 $\widehat{y}$ 表示 $Y$ 出现或未出现在其窗口中的父母状态。 我们还有 LR 测试仪：

(12) $$\text{LR score}=F_{{\chi }_{2^{|Z|}}^2}(-2[log{L}^{\ast }(X|Y,Z)-log{L}^{\ast }(X|Z)])$$

这里，$F(\cdot )$是具有$2^{|Z|}$自由度的卡方随机变量的累积分布函数。 然后，我们为每个测试者应用一个阈值$\tau$，即当分数小于$\tau$时，类型$X$与类型之间不存在因果关系$Y$。

4.3.2. 基于分数的方法

同样，有一些基于分数的方法可以应用于 GEM。 Bhattacharjya 等人 (Bhattacharjya 等人, 2018) 提出了 PGEM - 一种假设其强度函数仅受父类型是否发生在最近某个时间窗口内的影响的模型。 此外，他们使用条件强度${\lambda }_{x|u}$上的BIC标准来搜索每种事件类型的最佳父集，即推断其PGEM模型中的图结构。 图结构是不同类型事件之间因果关系的表示。

4.3.3. 传递熵

回想一下，传递熵 (TE) 可用于发现离散时间情况下的因果关系。 在这里，我们还可以将 TE 应用于事件序列（即点过程）来识别连续时间因果关系。 具体来说，Spinney 等人 (Spinney 等人, 2017) 构造了一个连续时间成对的转移熵：

(13)

其中$N_X$是目标进程中的事件数量，$\tau$是保持相应强度函数和 当 ${𝐓}_{Y\to X}=0$ 时，进程是独立的。 我们可以类似地定义条件TE。 现有一些一致的方法来估计连续时间 TE 及其条件形式（Shorten 等人，2021）。

5. 应用领域

时间因果发现已广泛应用于许多领域，例如科学事业（地球科学（Runge等人，2019a）、神经科学（Reid等人，2019；Weichwald和Peters，2021； Siddiqi 等人, 2022)，生物信息学(Sachs 等人, 2005))，工业实现（异常检测(Qiu 等人, 2012)，根本原因分析（Vuković 和 Thalmann，2022；Liu 等人，2021；Assaad 等人，2023）、在线系统中的商业智能（Arabzadeh 等人，2018）、视频分析(Yi 等人, 2020))。 表5总结了应用领域和相应的研究。 对于科学研究来说，习得的因果关系通常不应被视为最终结果，而应被视为进一步研究的起点和假设（Mäkelä 等人，2022）。 作为促进者，因果发现可以在工业环境中的多阶段方法中发挥支持作用（Vuković 和 Thalmann，2022）。 在本节的其余部分中，我们将回顾包括地球科学、异常检测和根本原因在内的三个领域，分别解释将时间因果发现纳入科学努力和工业实施的主要工作流程。

表 5. 主要研究时间因果发现应用。

Groups	Application areas	Studies
Scientific endeavors	Earth science	Climate change detection and attribution (e.g., (Lozano et al., 2009b)); Quantifying climate interactions (e.g., (Runge et al., 2014)); Latent driving force detection (e.g., (Trifunov et al., 2019; Shadaydeh et al., 2019)); Causality validation between temperature and greenhouse gases (e.g., (Van Nes et al., 2015)).
	Neuroscience	Dynamic causal models for neural connectivity (e.g., (Penny et al., 2004, 2010; Jafarian et al., 2020)); Granger causal models for neural connectivity (e.g., (Kaminski et al., 2001; Stokes and Purdon, 2017; Sheikhattar et al., 2018; Kim et al., 2011)); Causal inference from noninvasive brain stimulation (e.g., (Bergmann and Hartwigsen, 2021)).
	Bioinformatics	Modeling gene regulatory network (e.g., (Li and Ngom, 2013; Li et al., 2016; Verny et al., 2017; Patil and Vaida, 2022; Wu et al., 2022a)).
Industrial implementations	Anomaly detection	Causal structure as detection reference (e.g., (Qiu et al., 2012; Behzadi et al., 2017; Apte et al., 2021; Yang et al., 2022a)); Detection from imbalanced data (e.g., (Huang et al., 2020c)).
	Root cause analysis	Oscillation propagation tracing in the control loop (e.g., (Landman et al., 2014; Landman and Jämsä-Jounela, 2016; Chen et al., 2017; Lindner et al., 2018)); Alarm flood reduction (e.g., (Wang et al., 2015; Rodrigo et al., 2016; Wunderlich and Niggemann, 2017)); Industrial knowledge combined analysis (e.g., (Landman and Jämsä-Jounela, 2016; Cao et al., 2022; Thambirajah et al., 2009; Winchester et al., 2022)).
	Business intelligence in online systems	User interest prediction (e.g., (Arabzadeh et al., 2018; Hauffa et al., 2019)); Social media analysis (e.g., (Chang et al., 2013; Tsapeli et al., 2017; Kuzma et al., 2021; Chen et al., 2020); Online advertising (e.g., (Nuara et al., 2019; Yao et al., 2022b; Chu et al., 2020)); User-item interaction in recommendation (e.g., (Shang and Sun, 2020)); User activity modeling (e.g., (Li et al., 2017; Yao, 2022)).
	Video analysis	Video analysis and reasoning (e.g., (Yi et al., 2020; Li et al., 2020a)); Interpretable Gait Recognition (e.g., (Balazia et al., 2022)).
	Urban data analysis	Trajectory pattern mining (e.g., (Chu et al., 2016; Yang et al., 2022b)); Traffic flow prediction (e.g., (Li et al., 2015)); Visual urban and causal analytics (e.g., (Deng et al., 2022a)).
	Clinical data analysis	Causal chain discovery (e.g., (Wei et al., 2022)); Hypothesis testing (e.g., (Pandey, 2021)); Stable causal structure learning (e.g., (Rahmadi, 2019)).
	Signal processing	Blind source separation (e.g., (Testi et al., 2020; Testi and Giorgetti, 2021)); Compressed sensing (e.g., (Kathpalia and Nagaraj, 2022)).
	Financial analysis	Causal discovery for financial news (e.g., (Tetereva, 2018; Rambaldi et al., 2015)).
	Military	Battlefield sequential events analysis (e.g., (Li et al., 2022)).
	Robotics and dynamic control systems	Identifying causal structure (e.g., (Baumann et al., 2020)); Causal generalization (e.g., (Sheikhlar et al., 2021)).

地球科学和气候变化研究： 时间因果发现方法已广泛应用于地球科学和气候变化研究领域(Lozano 等人, 2009b; Ebert-Uphoff and Deng, 2012; Runge 等人, 2014; Van Nes 等人, 2015; Hannart等人，2016；龙格等人，2019a；特里福诺夫等人，2019）。 气候是一个复杂而混沌的系统，包含时空信息。 由于不确定性、简化以及与观测数据的差异，基于正演模拟的传统气候模型在描述此类系统时存在固有的局限性（Lozano等人，2009b）。 然而，常用的以数据为中心的方法，如滞后互相关和回归分析，旨在深入了解气候过程之间的相互作用机制，可能会导致该领域得出模糊的结论（Runge等人，2014） 。 为了克服上述问题，有意义地描述感兴趣的参数之间的因果关系并做出断言是合理的。 具体来说，（Lozano等人，2009b）提出了基于群弹性网络的时空格兰杰建模来进行气候变化检测和归因，其中极值理论对气候中的极端事件进行建模和归因，例如严重的热浪和洪水。 在（Runge等人，2014）中，提出了一种图形格兰杰模型和因果相互作用强度测量，以量化气候相互作用的强度和延迟，并克服普通相关或回归方法可能出现的伪影。 另一个挑战是存在未观察到的混杂因素，这可能导致错误的归因或成为不可忽视的驱动因素。 一系列工作（Trifunov等人，2019；Shadaydeh等人，2019）通过估计混杂变量之间的因果联系强度来检测气候异常事件的潜在驱动力。 此外，在气候系统中，一些感兴趣的参数表现出强耦合性，从而给因果方向的识别带来了困难。 (Van Nes 等人, 2015)采用专为强耦合动态系统设计的收敛交叉映射（CCM）技术来识别温度与温室气体之间的因果关系，其中统计关联性有据可查，但从观察到的数据中提取的因果关系却不同。 (Runge 等人, 2019a) 还提供了地球系统中时间序列因果发现的最新概述，其中描述了方法开发和科学努力的未来工作途径。

工业时间异常检测： 在工业系统中，检测来自传感器、日志、物理测量、系统设置等的海量时间数据中的异常是有意义的，同时也具有挑战性。 这些异常可以大致分为单变量异常（已被广泛研究）和依赖性异常（检测起来更具挑战性，但在现实应用中很常见）。 由于挑战主要来自数据的高维度和复杂依赖性，方法(Qiu 等人, 2012; Behzadi 等人, 2017; Apte 等人, 2021; Yang 等人, 2022a; Huang 等人, 2020c) 基于时间因果发现，通过提供高效、稳健和可解释的结果，在依赖性异常检测中发挥了不可忽视的作用。 因果发现可以促进底层系统生成机制的检测。 该系列工作的关键思想是首先从多元时间序列构建因果图，然后根据提取的因果关系检测异常。 具体来说，(Qiu 等人，2012；Behzadi 等人，2017)中，分别在参考数据集和测试数据集上建立了 Granger 图模型，分布差异（如 KL两个学习模型之间的散度和 Jensen-Shannon 散度）被计算为异常度量。 在(Apte等人, 2021)中，如果基于格兰杰因果关系的推断关系违反了领域知识或经常观察到的形式，则称为因果异常。 最近，（Yang 等人，2022a）也采用因果视角来检测多元时间序列异常并在 AIOps 应用中利用。 在这项工作中，计算成本降低了，因为它不是直接对联合分布进行建模，而是根据学习的因果结构对分解的分布模块进行建模，其中每个模块对应于局部因果机制。 此外，针对异常数据点较少的不平衡飞行数据，（黄等人，2020c）提出了一种基于非线性格兰杰因果学习的时间序列分类方法。

制造过程根本原因分析： 根本原因分析是确保工业环境中过程安全和生产力的一项重要任务，其中制造过程是临时且复杂的场景，通常由多个过程单元和大量反馈控制回路组成。 然而，由于对公平性、问责制和透明度（又称 FAT 原则（Shin 和 Park，2019）)的要求不断提高，该领域强大的机器学习方法的接受受到阻碍，特别是在敏感用途中案例（Vuković 和 Thalmann，2022）。 为了缓解这个问题，提取因果关系等知识在该领域至关重要。 过去十年见证了根本原因分析的因果发现方法的激增（Landman 等人，2014；Wang 等人，2015；Rashidi 等人，2018；Liu 等人，2021；Vuković 和 Thalmann，2022） 。 例如，利用格兰杰因果关系、转移熵及其变体等时间因果发现方法来追踪控制环中的振荡传播（Landman 等人，2014；Landman 和 Jämsä-Jounela，2016；Chen 等人，2017；林德纳等人，2018）。 警报泛滥被认为是工业事故的主要原因，减少警报泛滥是工业根本原因分析的另一个方面。 在三种典型的滋扰警报（即重复警报、长期警报和后果警报（Henningsen and Kemmerer，1995）)中，抑制后果警报并在异常发生的情况下提供适当的警报是一项挑战。发生并传播。 识别警报之间的所有因果关系是有帮助的（Hollender and Beuthel，2007），并且是一系列工作（Wang 等人，2015；Rodrigo 等人，2016；Wunderlich 和 Niggemann， 2017） 在这项任务中利用因果发现方法。 此外，信息流、能量流等深厚的工业知识可以与因果发现相结合，消除虚假关系(Landman and Jämsä-Jounela, 2016; Cao 等人, 2022; Thambirajah 等人, 2009).

6. 绩效评估

在本节中，我们概述了时间因果发现中使用的基准数据集和评估指标。

6.1. 数据集

我们将简要介绍一些用于时间因果发现的数据集，包括 MTS 数据集和事件序列数据集。

MTS 因果发现的数据集范围从健康数据到财务数据。 我们讨论一些常用的数据集，这些数据集是公开可用的，并且具有因果图的基本事实。

• Lorenz-96模拟数据：它是（Lorenz，1996）中制定的用于模拟气候动态的非线性模型。 $d$ 维洛伦兹模型中的连续动力学由 $\frac{\partial {𝐱}_i^t}{\partial t}=-{𝐱}_{i-1}^t({𝐱}_{i-2}^t-{𝐱}_{i+1}^t)-{𝐱}_i^t+F,i\le i\le d$ 给出。 对于较高的强迫常数$F$值，系统动力学变得越来越混乱。作为标准基准，它被 使用（Tank 等人，2022；Khanna 和 Tan，2020；Marcinkevics 和 Vogt，2021；Chu 等人，2020；Li 等人，2023）。

• 线性VAR模拟数据：根据线性VAR模型生成时间序列测量值。 在（Tank等人，2022；Khanna和Tan，2020）中，它用于分析当真正的底层动态是线性时方法的性能。

• CMU人体动作捕捉（CMU MoCap）数据：来自CMU MoCap数据库的数据集，包含有关关节角度、身体位置的数据。 可以利用因果发现方法来提取身体不同区域之间的非线性依赖性（Tank等人，2022）。

• DREAM-3 in Silico Network Inference Challenge：在DREAM-3 IN Silico Network Challenge（Prill等人，2010）中，使用连续基因表达和时间序列数据来模拟时间序列数据监管动态。 从记录的基因表达水平轨迹推断出五个基因调控网络。 该数据集已用于评估（Tank等人，2022；Khanna and Tan，2020）中的因果发现算法。

• 血氧水平依赖性 (BOLD) 成像数据：在此数据集中 (Smith 等人, 2011)，BOLD 信号的时间排序样本测量人类受试者感兴趣的不同大脑区域。 它是使用动态因果建模功能磁共振成像 (fMRI) 正演模型生成的。 在（Khanna and Tan，2020；Nauta等人，2019）中，应用因果发现方法基于BOLD成像数据来估计人脑中的连接。

• 模拟金融时间序列：数据集 (Kleinberg, 2013) 使用因子模型创建，根据三个因子和一个投资组合来描述投资组合的回报-特定错误术语。 因此，真实的关系是已知的。 它由 (Nauta 等人, 2019) 使用。

至于事件序列，数据集范围从在线行为到电力。 然而，并非在所有情况下都能获得有关因果关系的真实信息。

• MemeTracker：这是一个数据集，捕获在线文章的网站、发布时间以及其中的所有超链接。 该数据集最初代表了模因在不同网站上的流动方式。 网站的域名和发布时间被视为事件类型及其发生时间。 而不同网站之间的超链接可以被视为因果关系的基本事实。 它由 （Achab 等人，2017；de Figueiredo 等人，2018；Zhang 等人，2020） 使用。

• IPTV观看记录：该数据集（罗等人，2014）记录了用户的观看行为，即在IPTV系统中观看的节目和时间。 节目类型和观看节目的时间可以分别被视为事件类型和事件发生时间。 它由（Xu等人，2016；Chen等人，2022a；Zhang等人，2020）使用。 但是，该数据集中不包含真实的因果关系。

• 电网故障事件数据：该数据集包括相量测量单元（PMU）内电压或电流信号的突变以及每个 PMU 的 ID。 使用该数据集进行因果诊断任务的任务是推断网格内的因果关系（Idé等人，2021）。 由于出于隐私考虑而未给出网络拓扑，因此这是一项非真实任务。

• 七国集团债券：该数据集（Demirer等人，2018）包含七国集团国家主权债券的每日回报波动率。 处理此数据集的目标是发现主权债券下面的因果网络（Jalaldoust 等人，2022）。 来自该领域的专业知识可以被视为基本事实。

6.2. 评估指标

在这一部分中，我们将解释文献中使用的不同指标。 给定以 $thre\in (0,1)$ 为阈值的边 $p(A_{ij})$ 的推断概率、因果图 中的地面实况边集以及因果图 $E_{GT}=\{(i,j):A_{ij}^{\ast }=1\}$ 中的地面实况缺失边集因果图$E_{MS}=\{(i,j):A_{ij}^{\ast }=0\}$，常用指标的定义和描述如下：

• 真阳性率（TPR）：作为因果发现结果中发现的公共边和地面实况邻接点与地面实况边总数的比率，TPR 指标定义为 $TPR=\frac{|\{(i,j):p(A_{ij})\ge thre\}\cap E_{GT}|}{|E_{GT}|}$.

• False Positive Rate (FPR)：与TPR类似，FPR是指因果发现结果中发现的共同边与groundtruth缺失邻接点占groundtruth缺失边数量的比例，即定义为$FPR=\frac{|\{(i,j):p(A_{ij})\ge thre\}\cap E_{MS}|}{|E_{MS}|}$。

• 接收器运算曲线下面积 (AUROC)：接收器运算曲线 (ROC) 的定义是，阈值 $thre$ 在 $0$ 和 $1$ 之间变化时，TPR 和 FPR 的比率。 ROC 下的区域 (AUROC) 随后被广泛用于评估因果发现算法的性能。

• 结构汉明距离（SHD）：SHD是一种度量，描述将发现的图转换为其对应的真实图所需的边版本数量，它总结了缺失边、额外边的数量和不正确的边缘。

7. 讨论和新观点

在本节中，我们首先讨论挑战和实际考虑因素，包括非平稳性、异质性、未观察到的混杂因素、二次抽样和专家知识。 然后，提供了时间因果发现的两个新视角，我们认为这将是未来研究的一个有希望的途径。

7.1. 挑战和实际考虑

数据的非平稳性： 在实际场景中，我们经常面临非平稳性，时态变量的概率分布以其原因为条件，甚至因果关系可能会随着时间的推移而变化，特别是对于时态数据。 在这种情况下，假设固定因果模型的因果发现方法可能会给出误导性的结果。 然而，多种类型的研究表明，非平稳性包含了因果发现的信息（Tian and Pearl，2001；Peters 等人，2016；Zhang 等人，2017；Huang 等人，2019）。 因此，正确解决应用中的非平稳性非常重要。 非平稳性可能是由底层系统的变化引起的，可以看作是大自然的软干预（Korb等人，2004）。 遵循这一想法，一系列工作（Zhang 等人，2017；Huang 等人，2020d）利用时间和域索引等代理来解释因果关系发生变化的非平稳性，并且提出了CD-NOD框架。 另一组研究没有利用信息非平稳性进行因果结构学习，而是专注于对时变关系进行建模（Gao 和 Yang，2022）。 此外，（Du和Xiang，2020）提出了缓慢变化的非平稳过程（例如演化谱和局部平稳过程）的方法。

数据的异质性： 在实际应用的因果发现中，数据的异质性体现在两个层面：（1）相互作用的时间过程是异质的（具有不同的分布），例如，不同站点的因果相关的气象观测分别受到几个主要天气系统的影响（Behzadi 等人，2019）。 （2）底层生成过程会随着数据集或不同领域的变化（Glymour等人，2019），例如不同市场的股票价格（Huang等人，2020d）或不同范式下的个体行为（Chu等人，2020）。 对于时间变量之间存在异质性的第一个条件，针对特定同质数据类型设计的传统因果发现方法的推断关系可能不准确。 作为补救措施，中提出了基于广义线性模型和最小消息长度等方法的格兰杰因果关系的几种变体（Behzadi 等人, 2017, 2019; Hlavácková-Schindler and Plant, 2020） 。 对于第二个条件，工作线(Zhang 等人, 2017; Huang 等人, 2020d)利用异质性的分布转变作为软干预来辅助因果结构学习，这与在非平稳数据中。 而第二种情况下的另一条因果发现方法（Chu 等人，2020；Löwe 等人，2022）侧重于在端到端框架内对异构数据中的典型结构进行归纳建模。

未观察到的混杂因素： 在实践中，我们经常会遇到违反因果充分性的情况，即，存在未观察到的混杂因素。 这种具有挑战性的设置可能会导致不正确的因果关系（Geiger等人，2015）。 正如表 3 所总结的，大多数时间因果发现方法无法以直接的方式处理未观察到的混杂因素。 设计了几种基于约束的方法，但没有因果充分性和方法。此外，通过在（Löwe等人，2022）中应用结构偏差来对未观察到的混杂因素进行建模。 最近几项被称为因果表征学习的研究对未观察到的混杂因素采取了新的视角。 将在(7.2.2)小节中详细介绍。

二次采样：在实际应用中，由于数据收集的困难，时态数据（尤其是时间序列）的采样率可能低于底层因果过程的速率。 用于子采样时间序列的普通因果发现算法可能会导致虚假因果关系和遗漏因果关系。 提出了几点评论和方法（Danks和Plis，2013；Gong等人，2015；Plis等人，2015；Gong等人，2017；Hyttinen等人，2016；Tank等人，2019）对于这个问题。

专家知识： 专业知识可以帮助实践中的因果发现过程。 融合专家知识的方法可以分为三种类型（Kitson等人，2021）：（1）软约束：学习过程会受到知识的影响(O'Donnell 等人, 2006)。 (2) 硬约束：学习到的结构必须符合强制要求(即，以概率$p=0$或$p=1$ 在（Asvatourian等人，2020）中，在具有时间依赖性暴露的结构学习中利用了硬约束。 (Sun 等人, 2021) 中的研究添加了禁止切片内依赖存在的先验知识，这有助于恢复未由先验知识显式编码的边缘。 （3）交互式学习：在学习过程中利用人类输入（Messaoud 等人，2009；Melkas 等人，2021；Zhu 等人，2022；Jin 等人，2021） 。

7.2. 新视角

7.2.1. 摊销和监督范式的扩展

在传统范式中，因果发现方法大多要么单独处理观测数据，要么为每个人训练一个不同的模型。 这些方法没有充分利用不同样本之间的共同结构，也没有充分利用因果结构清晰的数据集中的监督信息，因此存在小样本挑战和缺乏归纳能力等问题。 最近，在新的范式中进行因果发现来解决这个问题。 我们可以将它们大致分为两类：基于摊销建模（Chu等人，2020；Löwe等人，2022）的方法和基于监督的方法学习 （Benozzo 等人，2017；Wang 和 Kording，2022）。 我们在本小节中介绍它们，我们相信这是未来研究的一个有前途的途径。

在摊销模型中，为具有不同因果结构的个体训练全局因果发现框架。 对于时态数据的场景，这些方法已在 3.4.4 中详细介绍，作为格兰杰因果关系与归纳建模的深度学习扩展。 InGRA (Chu 等人, 2020) 利用原型学习来提取常见的因果结构，而 ACD (Löwe 等人, 2022) 提出了一种编码器-解码器框架来进行摊销因果发现。 这些方法充分利用海量样本的信息，能够推断新到达个体的因果关系，这在电子商务、社交网络和神经图像等现实世界的应用中非常有用。

另一条工作主要集中在将推理过程视为黑匣子，并通过监督学习学习从样本数据到因果图结构的映射。 这里的标签信息是因果结构，可以在合成数据集中轻松访问。 通过监督学习学习因果关系的早期工作（Lopez-Paz 等人，2015；Ton 等人，2021）仅限于学习成对因果方向，其中问题被转化为分类任务来区分$X\to Y$ 和 $Y\to X$ 使用观察到的样本。 后来在（Li等人，2020b；Petersen等人，2022）中扩展到发现图结构。 由于训练的标记信息通常源自已探索的合成数据或现实世界数据集，因此不能保证训练和测试数据的分布匹配或高度重叠的监督方法的要求。 在（Ma等人，2022；Ke等人，2022）中，利用邻位图和元学习等方法进行监督因果发现来解决这种“领域转移”问题。 对于时间设置，（Benozzo 等人，2017）提出了时间序列之间格兰杰因果关系的监督估计。 作为最近的进展，（Wang 和 Kording，2022） 中提出了一种学习因果发现的方法，其中在时间数据集上进行测试时，从具有已知因果关系的大型数据集中学习的算法优于传统范式中的算法例如功能磁共振成像。

7.2.2. 因果表示学习的扩展

从深度学习黑匣子中显式或隐式地提取特定现象的原因可能有利于下游任务。 上述因果发现方法侧重于推断观察变量之间的关系，或者从预先给出因果变量的前提开始。 尽管有些方法可以在未观察到的变量下学习因果关系。 存在现实世界的观察（例如，传感器测量、视频中的图像像素），这些观察一开始就没有很好地构造因果变量。 作为从观测变量中发现因果关系的一般化，因果表征学习（Locatello 等人，2020 年；Schölkopf 等人，2021 年；Yang 等人，2021 年）最近越来越受到关注，其目的是学习底层系统中因果因素的表征。 它根据观察估计潜在因果变量图。

因果表示学习领域的一系列工作基于解缠结和独立成分分析（ICA）来识别变化的独立因素。 该方法的核心是假设相互独立的潜在因素。 识别真正的潜在变量是很困难的，尤其是在一般的非线性情况下。 作为一种补救措施，最近的方法（Locatello等人，2020；Khemakhem等人，2020；Hyvärinen和Morioka，2017，2016）利用多个视图、辅助变量或时间结构中的附加信息，并结合深度学习方法，例如 VAE 和对比学习。 最近在 （Gresele 等人，2021；Monti 等人，2019）中建立了 ICA 和因果关系之间的联系。 在时间数据的背景下，在潜在时间因果过程估计(LEAP)（Yao等人，2022c）中讨论了时间序列因果变量的可识别性。 它首先提供非参数、非平稳设置和参数设置中的因果可识别性条件。 然后，它提出了一个提取潜在因果关系的学习框架，通过强制执行假设条件，用学习的因果过程网络扩展了 VAE。 由基于流的估计器建模的非平稳噪声可以被视为帮助识别的软干预。 与 LEAP 一致，后续工作（Yao 等人，2022a）将识别理论扩展到更一般的情况。

另一项工作利用干预和数据增强来帮助识别潜在的因果关系。 在数据增强下，（von Kügelgen等人，2021）证明了常见的对比学习方法可以分块识别保持不变的因果变量。 对于时间设置，提出了CITRIS(Lippe等人，2022b)。 这是一个学习因果表示的 VAE 框架，其中潜在的因果因素可能已被干预。 通过使用干预目标信息进行识别，CITRIS 不受功能或分布形式的限制。 此外，CITRIS 中的因果因素被视为标量或潜在的多维向量，这在复杂场景中更实用。 沿着这一工作思路，在 iCITRIS （Lippe 等人，2022a） 中提取瞬时因果关系。

8. 结论

时态数据中的因果发现对于理解动态和估计感兴趣的因果效应至关重要。 本文回顾了时间因果发现的两类：多元时间序列因果发现和事件序列因果发现。 多元时间序列因果发现可以分为四类，包括基于约束、基于分数、基于 FCM 和 Granger 因果模型。 回顾了每种类型的主要思想和最新进展。 对于事件序列中的因果发现，我们可以将这些算法分为基于约束、基于评分和格兰杰因果模型，这些模型与多元时间序列因果发现相一致。 我们注意到，由于格兰杰因果关系和霍克斯过程之间的自然匹配，格兰杰因果模型特别适合事件序列。 为了弥合丰富的时间因果发现算法与现实世界影响之间的差距，我们介绍了几项主要研究，包括科学努力和工业实施。 我们还提供了广泛的资源列表，包括数据集和指标，可用作该领域未来研究的指南。 虽然许多算法都提供了理论或经验保证，但推断关系的质量取决于许多问题，包括非平稳性、异质性、未观察到的混杂因素、二次采样和专家知识。 我们讨论这些挑战和实际考虑因素。 最后，我们介绍了因果发现的新视角，其中描述了摊销建模、监督学习和因果表示学习方面未来工作的途径。