Identifying Causal Structure in Dynamical Systems
Authors: Dominik Baumann
Friedrich Solowjow
Karl Henrik Johansson
Sebastian Trimpe
Venue: Accepted by TMLR ()

Abstract: Mathematical models are fundamental building blocks in the design of dynamical control systems. As control systems are becoming increasingly complex and networked, approaches for obtaining such models based on first principles reach their limits. Data-driven methods provide an alternative. However, without structural knowledge, these methods are prone to finding spurious correlations in the training data, which can hamper generalization capabilities of the obtained models. This can significantly lower control and prediction performance when the system is exposed to unknown situations. A preceding causal identification can prevent this pitfall. In this paper, we propose a method that identifies the causal structure of control systems. We design experiments based on the concept of controllability, which provides a systematic way to compute input trajectories that steer the system to specific regions in its state space. We then analyze the resulting data leveraging powerful techniques from causal inference and extend them to control systems. Further, we derive conditions that guarantee the discovery of the true causal structure of the system. Experiments on a robot arm demonstrate reliable causal identification from real-world data and enhanced generalization capabilities.

Published in Transations on Machine Learning Research (07/202)
Identifying Causal Structure in Dynamical Systems
Dominik Baumann
dommnnkk. boumann Oit.
Deportment of Information Technology
Uppsala University
Uppsala, Sweden
Friedrich Solowjow
Institute for Data Science in Mechanical Engineering RWTH Aochen University
Aochen, Germany
Karl H. Johansson
kellej@kth.s
Division of Decision and Control Systems and Digital Futures KTH Royal Institute of Technology
Stockhoim, Sweden
Sebastian Trimpe
trimppe@dsmeerrtt-acckee
Institute for Data Science in Mechanical Engineering RWTH Aochen University Aachen, Germany
Reviewed on OpenReview: https: // openreview. net/ forum? id= X2BodlyIof
Abstract
Mathematical models fundamental building blocks in the design of dynamical control systems. As control systems are becoming increasingly complex and networked, approaches for obtaining such models based on first principles roach their limits. Data-driven methods provide an alternative. However, without structural knowledge, these methods are prone to finding spurious correlations in the training data, which can hamper generalization capabili- ties the obtained models. This can significantly lower control and prediction performance when the system exposod unknown situations. A preceding causal identification can prevent this pit.fall. In this paper, we propose method that identifies the causal ructure of control systems. We design experiments based on the concept of controllability, which provides a systematic way to compute input trajectories that steer the system to specific regions in its state space. We then analyze the resulting data leveraging powerful techniques from causal inference and extend them to control systems. Further, we derive conditions that guarantee the discovery of the true causal structure of the system. Experiments on a robot arm demonstrate reliable causal identification from real-world data and euhanced generalization capabilities.
1 Introduction
When learning models for dynamical control systems, we ideally would like to obtain models that (i) gener- alize well to new domains, (ii) are interpretable, and (iū) computationally efficient. However, generalization to new domains and interpretability are particular weaknesses of current black-box machine learning meth- ods (Schōlkopf, 2022 Rudin 2019). We address this problem by first. identifying control system's causal structure and subsequently using this structural knowledge for model learning.
As the understanding of causality differs depending on the domain, we first provide some intuition of what kind of structure seek to identify. For this, following Fichler (2022), consider two notions. The first is

---

Published in Transations on Machine Learning Research (07/202)
temporal precedence: causes precede their effects. Temporal precedence is the typical causality notion used in systems theory (Hespanha 2018 p. 31). However, here, we mainly focus on the second notion, physical inftuence: manipulating causes changes the effects. In other words, we seek experimental routines and tests that enable control systems to learn (i) what is the influence of their internal states on one another, and (it)
which of their inputs influence which internal stales.
Incorporating the causal structure into the model learning problem targets the three characteristics of ideal models outlined above. A key reason for the shortcomings of most current machine learning algorithms concerning those characteristics is that they mostly rely pure statistical analysis of the data (Pearl 2018). Thus, when considering stochastic systems and finite sample sizes, these algorithms will, due to spu- rious correlations, typically find connections between variables, although those variables do not influence one another based on the underlying physical equations. This can lead to catastrophic errors when extrapolating outside the training data, thus, diminishing generalization capabilities. Further, such models are less inter- protable since is unclear whether connections between variables in the model indicate a causal influence or merely spurious correlation found in the training data. A causal analysis mitigates both problems since it identifies which variables actually have causal influence each other. Incorporating this knowledge into the model learning problem yiclds more interpretable models that generalize better to new domains (Pearl & Mackenzie 2018 Lastly, when only connections between variables that actually have a causal influence on one another are considered, this also reduces the model's parameter space. Therefore, the causal analysis also helps make the model computationally more efficient.
In this paper, we propose an algorithm that automatically identifies the causal structure of a control system. The problem of inferring causal structure from observational data, i.c., given data for which we cannot influence the data generating process, has been addressed using various methods, see [Spirtes et al. (2000); Peters et al. (2017) for an overview. For control systems, solely relying on observational data is not necessary. Such systems are equipped with an input, which they can use to actively conduct experiments for causal inference. Causal inference from experiments, or interventions, has also been studied, most prominently in the context of the do-calculus [Pearl] 1995). However, there it is assumed that state variables can be directly influenced by the input, which is often not possible in control systems. In control systems, it is essential to consider a proper notion of controllability, i.c., how the system can be stecred to particular regions in the state space through appropriate input trajectories. To the best of our knowledge, such notions have not been considered yet. in the causality literature. In the literature on model learning and system identification, we find many techniques for learning mathematical models for control systems (Ljung, 1999 N Suyyen-T Tuong& & Peters 2011 Schoukens & Ljung. 2019). In system identification, the problem of the exploding number of parameters for black-box methods has, for instance, been addressed using regularization techniques (Schoukens & Ljung. 2019). While these methods can roduce the number of parameters, they may exclude parameters representing causal influence or include parameters representing spurious correlation depending on the regularization parameter. Thus, the obtained models may fall short generalization and interpretation capabilities.
Contributions. In summary, existing methods from causal inference do not consider a proper notion of controllability, while existing methods from system identification cannot give formal guarantees on finding the true causal structure. In this paper, we bridge this gap and present an algorithm that identifies control system's 's causal structure through an experimental design based on a suitable controllability notion and a subsequent causal analysis, for which can provide formal guarantecs. For the causal analysis, we leverage powerful kernel-based statistical tests based on the maximum mean discrepancy (MMD) (Gretton et al. 2012). Since the MMD has been developed for independent and identically distributed (f..d.) data, we extend it by deriving conditions under which the MMD still yields valid results and by coming up with a test statistic for hypothesis testing, despite non-i.id. data. In terms of controllability, we investigate three different settings: (i) exact controllability, where we can exactly steer the system to a desired position, (wi stochastic controllability, where we can only steer the system to an e-region around the desired position, and (iii the special case of linear systems with Gaussian noise that are controllable in the sense of Kalman [Kalman 1960a), as they represent a widely studied class of systems. We demonstrate the proposed method's applicability by automatically identifying the causal structure of a robotic system and a simulated, nonlin- ear quadruple tank system. Further, we show improved generalization capabilities for the robotic system
2

---

Published in Transations on Machine Learning Research (07/202)
and reduced computational complexity for the quadruple tank system, both inherited through the causal identification.
2
Related Work
To the best of our knowledge, no other algorithm seeks to identify the causal structure of dynamical control system through experiments based on a suitable controllability notion. However, several works in causal inference aim to infer the causal structure of dynamical systems, and several methods in system identification seek to reduce the parameter space or identify structural properties of control systems. In this soction, we discuss those works.
Causal inference for dynamical systems. Causal inference in dynamical systems or time series has been
studied in Demiralp & Hoover (2003); | Eichler (2010); Monela et al. 2001Tp; Malinsky & Spirtes (0018) using
vector autoregression, in Peters et al. 2013) based on structural equation models, in
Entner & Hoyer (2010).
using the fast causal inference algorithm (Spirtes et al. 2000), and in (alvi et al.
(202]) and Quinn et al.
(2011), applying kernel mean embeddings and directed information, respectively. In Pfister et al. (2019): the
authors develop a procedure for learning the causal structure of kinetic systems. A more extensive overview of causal inference methods that can be applied to time-series data, with a particular focus on Earth system sciences, is provided in Runge al.] (2019a). While all these methods make causal inference for dynamical systems, none of them investigates experimental design. Instead, they aim to infer the causal structure from obscrvational data and, thus, need additional assumptions to arrive at statements about the causal structure. Dynamical systems, as considered in this work, can be actively influenced through a control input. Hence, we can design experiments and do not need to rely on the data being sufficiently rich.
Experimental design. A well-known concept for causal inference from experiments is the do-calculus. In the basic setting, a variable is clamped to fixed value, and the distribution of the other variables conditioned this intervention is studied (Pearl 1995 Extensions to more general classes of interventions exist, see, e.g., Yang et al. (2018) Shanmugam el al. 2015), but they consider static models, which is different from the dynamical systems studied herein. Causal inference in dynamical systems or time series with interventions has been investigated in Eichler (2012 Peters et al.] (2022). : Mooij et al. (2013): Rubenstein et al. (2018):
Sokol & Hansen (2014). However, therein is assumed that one can directly manipulate the variables, e.g.; by setting them to fixed values or forcing them to follow a trajectory, which is typically impossible
practice. None of those works considers various degrees of controllability, which is the case in this paper.
Thus, they are not readily applicable to control systems.
Model selection and regularization. As an alternative to directly testing causal relations between variables, several methods exist that identify a dynamic model, trading off model complexity and accuracy. Typically, this is done by letting the algorithm select from a set of candidate models. In system identification, the Akaike information criterion (Akaike 1973) and the Bayesian information criterion (Schwarz 1978) are two well-known examples of such methods. In neuroimaging, there are dynamic causal models Friston et al. & Lipson 2007| Schmidt & Lipson 2009 Brunton et al.) 2016). In all cases, the true causal structure of the system can only be revealed if a model representing this structure is part of the candidate models. Further, they typically use a regularization parameter to find a trade-off between model complexity and accuracy. This parameter punishes model complexity (c.g., the number of paramcters) and rewards goodness of fit. Thus, it also depends on the specific choice of this regularization parameter whether or not these algorithms find a model representing the system's true causal structure.
2003 Stephan et al., 2010). A third family of methods are symbolic regression techniques [Bongard
Structure detection in dynamical systems. Revealing causal relations in dynamical system can be interpreted as idcntifying its structure. Related idcas exist in the identification of hybrid and piccewise affine systems (Roll et al. 2004 Lauer & Bloch 2018). These approaches try to find a trade-off between model complexity and fit but cannot guarantee to find the true causal structure. Further methods that identify structural propertics of dynamical systems can be found in topology identification (Materassi Innocenti] 2010 Shahrampour & Preciado 2014 van den Hof et al., 2013) and complex dynamic networks (Boccaletti)
et al. 2006: Liu et al. 2009 Yu 2010). Those works seek to find interconnection between subsystems
3

---

Published in Transations on Machine Learning Research (07/202)
instead of identifying system's inner structure done herein. Moreover, while the mentioned works often rely on restrictive assumptions such as known interconnections or linear dynamics, our approach can deal with nonlinear systems and does not require prior knowledge.
Kernel mean embeddings. For causal inference, we will leverage concepts based on kernel mean embed- dings, which have been widely used for causal inference (Peters et al. 2017 Chen et al.) 2014 Lopez-Par et al. 2015). A downside of those methods is that they typically assume that data has been drawn i.i.d. from the underlying probability distributions. Extensions to non-ii.d. data exist (Chwialkowski & Gretton 2014 Chwialkowski et al. 2014), but rely on mixing time arguments. Dynamical systems, as investigated in this work, often have large mixing times or do not mix at all (Simchowitz et al. 2018). Therefore, these types of analyses are not sufficient in this case.
3 Problem Setting and Main Idea
We consider dynamical control systems of the form
(t)) = ((0)))
(1)
with discrete time index t E N, dynamics function f. state x(t) E X' C R", state space X, input u(t) e U C IR"", imput space U, and e(t) E R an independent random variable sequence. The notation (0:t) here denotes the sequence u(0), u(1)..... u(t). In the following, we will omit this notation and simply write u or v if we consider the whole trajectory. The description of the system in equation 1 is different from the standard, incremental version ztt + 1) f(x(t), u(t), v(t). Specifically, equation emphasizes the dependence of x at time f on the initial state. Equation can readily be obtained by itcrating f starting from r(0).
We will lend notation from the do-calculus Pearl 1995 for causal inference. In this notation, P(x(()) (o) =0)) defines the conditional probability distribution æx() given that we set the initial condition of the ith component of the state vector x to a specific value (0) while all other components of x and the inputs u are left unaffected. In this article, we show how such do-conditional distributions can be realized in dynamical systems when considering a suitable notion of controllability and how can use them for causal inference.
3.1 Problem Setting
The system description from equation is stochastic and, hence, induces a probability distribution P(x) over trajectories æ. Based on this, we define non-causality adapting standard notation from the do-calculus Pearl 1995) to dynamical systems as in equation [
Definition 1 (Global non-causality). The state variable zj does not cause x; Pf =0)
P(z, (0)zz((0)= = 2!! (0))) for all (0) and and (0). The superscripts and I1 denote two eaperimental designs, ie.. different choices of (0). Similarly, uj does not cause I; if P(x;; do(u; = w/)) = PP(x I do(u; = wff)) for all ut and w.l.
The main objective of this article is to develop an algorithm to test for nou-causslity in the sense of definition I with theoretical guarantees. That is, we seek to test whether the initial condition of x, respectively the input trajectory uj changes the probability distribution of the resulting x; trajectory. The do-operator denotes that we force ž, to specific initial condition and uy to be specific input trajectory, while all other initial conditions and input trajectories remain fixed between experiments. While this is a common assumption of the do-calculus, in dynamical systems, we camnot just set state variables to specific values. Instcad, we need to consider a proper notion of controllabilily. In the following, we discuss how the do-operator can be applied to dynamical systems, i.e., how we can, based on an appropriate controllability notion, steer state variables toward specific values as required test for non-causality. Subsequently, we develop an algorithm that identifies for each pair 2; and zy and for each pair x; and dy of control system as in equation weether not both variables are non-causal according to definition
4

---

Published in Transations on Machine Learning Research (07/2022)
(0)
Exp I
2
2) 1...m)y(() (0) = =x (0), wt w!ll
MMD
(0)
Exp II
Figure 1: Experimental design for causal inference. We design two experiments, where all initial conditions and input trajectories are constant except for z0))). If the resulting trajectories of xs are different in both experiments, we have evidence that the change in (0) caused this change.
3.2 Main Idea
For causal inference, we will exploit that we can influence equation [ through u. To make the notion of how we can influence the system precise, we adopt controllability definitions for stochastic systems Sunahara et al. 1974 Bashirov & Kerimov 1997/:
Definition 2. The system described by equation 7 is said to be completely e-controllable in probability rj in the normed square sense in the time interval [0, tel if for all desired states Ides and initial states z=0) from X, there erists an input sequence from U such that Pr[[[z(tr) - zdeell 2 c} S 1 77, with 0 < 1.
A variety of methods exist to identify or learn models for systems as in equation I e.g., Gaussian process regression (Williams & Rasmussen 2006) or fitting lincar state space models using least squares (Ljung 1999). In the following, we assume that we can obtain an estimate f of the investigated system in equation (including an estimate of the distribution of u(t) using techniques from system identification and model learning. The concrete method that is used to obtain f is independent of the developed causal identification procedure. However, in later sections, we will specify requirements the model est istte f needs to satisfy. This model estimate f. obtained without further assumptions or physical insights, will, due to the stochasticity of equation [] almost surely entail spurious correlation and suggest causal influences that are actually not present in the physical system. Nevertheless, it will allow us to (approximately) steer the system specific initial conditions and start experiments from there.
We propose two lypes of experiments to test for causal relations. For the first lype, we investigale whether zj causes Zj. We conduct two experiments (denoted by I and I1) with different initial conditions =²(() # (0). while all others are kept the same (cf. figure Dp. This can be formalized as
(0) = (0) for all ( 2)(0)) + (0), ui(() = (t) for all e, t.
(2a)
By comparing the resulting trajectories of 2! and x,", we can check whether the change in (0) caused a different behavior. For checking the similarity of trajectories, we will use the MMD, whose mathematical definition we provide in the next section. The second type experiments analogous to the first, but instead of varying initial conditions, we consider different input trajectories u/ + u
(2b) Note that for the testing experiments, we design open-loop trajectories. That is, the input cannot depend on the system's current state. This is essential since creates independent trajectories for which we can leverage the MMD as a similarily mensure.
(0) = (0) for all f, (t())= = (t) for all l + j. for all t, a)(() + (t) for all t.
Remark 1.
3.3 Road-map
We propose an algorithm that consists of two steps: (i) we design experiments, and (i) we analyze the resulting data using the MMD. The data obtained from these experiments must fulfill specific requirements allow for proper causal inference. In the next section, we provide the mathematical definition of the MMD and deduce the requirements the experimental data. We then subsequently develop a suitable
5

---

Published in Transations on Machine Learning Research (07/2022
experimental design based on the controllability notion stated in definition We state conditions under which this experimental design yields data from which the MMD can provably infer the true causal structure in the infinite sample limit and derive a hypothesis test for finitely many samples. Until here, we focus on a rigorous convergence analysis. However, this requires us to do many experiments on the real system. In section 5 we propose a heuristic algorithm that is more efficient in terms of the number of experiments and computations but forgoes some of the guarantees. Lastly, scction G we demonstrate the method's applicability on real robotic system and a quadruple tank process and present comparisons with sparse identification and causal discovery method on synthetic linear example.
4 Causal Identification for Dynamical Systems
We will now develop the causality testing procedure. First, we introduce the MMD (Gretton al. 2012 which we shall use as a similarity measure. The MMD can be used to check whether two probability distributions P and Q are equal based on samples drawn from these distributions. Let X and Y be samples drawn ii.d. from IP and Q respectively. Further, let H be a reproducing kernel Hilbert space (Sriperumbudur et al. 2010 with canonical feature map d: X + H. The MMD is defined as
(MD(P, Q) ||Ex~>
The feature map & can be expressed in terms of kernel function k(,;'). :): where k(x,y) = (0(x), ((y)) If the kernel is characteristic, we have MMD(P,, Q) 0 if, and only if, P = Q (Fukumizu et al. 2008 Sriperumbudur et al. 2011). In the remainder of the paper, we always assume a characteristic kernel (e.g..
the Gaussian kernel).
In general, also other measures that compore probability distributions may be used for our algorithm. An overview of such methods is provided in (Sriperumbuurr al.) (2012). popular crample would be the Kullback-Leibler divergence [Kulbbcck & Leibler 1951). In this paper, we propose to use the MMD since it allows us to compare probability distributions without actually estimating them, provides theoretical guarantees, and can be computed efficiently.
Remark 2.
In the following, we derive conditions that allow one to provably identify causal relations. We investigate three settings. First, we discuss the case where we can precisely steer the system to desired initial conditions (i.e., c = in definition 2 We then extend this to e > O, which requires a stricter controllability definition. Finally, we show that for linear systems with additive Gaussian noise, widely studied class of systems, the conditions stated by Kalman (Kalman 1960a) are sufficient, and the identification is substantially easier.
4.1 Exact Controllability
When considering control systems, instead of obtaining single i.i.d. samples from stationary distributions, we receive sequences of random variables sampled from a stochastic process as in equation [ This data is often non-i.i.d.. The objects of interest, whose distributions we want to compare, are then the z; trajectories obtained from two different experimental settings. To simplify notation, we denote trajectory obtained from the first setting by a! and the joint probability distribution over the trajectory states by pl. and equivalently for z! and plll. We sample from pd' and pill by designing two experiments as in equation 2 with fixed length T and repeating each experiment m limes. That is, we obtain 2m sequences of T random variables sampled at discrete intervals of fixed length (i.e,, t ++ always has the same length). These samples we denote by x]. x!! TR. Note that (i) all m sequenees are sampled from the same distributions Ipl and pll, (ii) while the distributions are non-stationary along time, the distributions for multiple experiments of fixed length T are identical, and (iii) all m sequences are independent of each other. The MMD then reads
MMD(x x = || (=x)))((((()) E [[([""" [d(x"))]]
(4)
If we now design experiments using equation 2 we can check the similarity of x trajectories with equa- tion 4 Following /zabó & Sriperumbudur (2018), , if we use characteristic kernel equation 4 the general properties of the MMD test are still valid given that the initial conditions x z²(0) and x z²(0), respectively the input trajectories ut and "ll, are i...d,, which is the case for our design. That is, MMD suggests that the
6

---

Published in Transations on Machine Learning Research (07/202)
distributions of the trajectories are different, we can conclude that there is causal influence. However, the other direction is less straightforward: for a system as in equation variables are dependent, as can be seen in the following example:
1 the MMD may be zero, even though
Example 1. Assume a control system with x (1 (1(t+)) = xn (f)az(t) and 12(t + 1) = u(t), and an input signal u(t) that is different from O. If we choose 21 (0) = 0, the trajectory of xı will, despite the fact that D2 clearly has a causal influence on 3ı, always be O no maller the inilial condilion v2(0).
To address this, we define the concept of local non-causality:
Definition 3 (Tocal non-causality).. Let Xnc C X and Ume U with Xuc U Une + B. The state variable 2j does locally not cause T, if P((,; do(z,(0) =)() = P((,, do(z)(0) = =0]) for all (0 and (0) given that the sequence x is entirely in Xne and the sequence u is entirely in Unn Similarly, the imput variable uj does locally not cause I; P(z; do(u; = uf)) = PP(;; do(u; = =uf)) for all u/ and uj given Il
that x is in Xue and u in Une-
The uon-causality becomes global if Xnc X' and Une u. To properly test for causal relations, we need to ensure that the experimental design in equation Siooonyyiee yields initial conditions and input trajectories that are not inside Xne and Une. For this, we propose to design experiments based the estimated model f. In particular, we utilize simulated trajectories based on the model f, which we denote by (& | f) and where for each t we have att)= f (a),,,, 0). We then need to make the following assumption about the system described by equation and the estimated model f:
i for which Xue Ule + Ø. Consider further two independent experimental designs following equation 2 with inátial conditions x0)),))1(0) and impul. trujectories w², nd for the first, and initial conditions x (0), L.V (0) and input trajectories w _ for the second experiment. Assume that all individual inputs u(t) that are part of u III or TV are in Une and that all states z(t) that are part of the simulatedd trajectories (E!!! f) or (œtV f) are inside Xuc. Further, assume that for the first experiment there erists some u(t) or some zz() that is not part of Une respectively Xnc. For such cases, we assume MMD(s!", a!V 1 < MMD(s), 2!! f)).
Assumption 1. Consider a dynamical system as in equation
This assumption does not require that f captures the causal structure. That is, if variable 2, does not cause variable T4, the model f may still suggest that they are causally related such modeling errors will then be accounted for through the identification procedure we propose in this paper. However, we require that f captures that the influence of r, on x; is lower in regions in which they are locally non-causal. Intuitively, the assumption says that. if at does not influence xy in certain parts of the state/ state/action space, our model may, due to spurious correlation, stiill assign some influence of I; on Iy in those regions. Nevertheless, we expect the model to assign stronger influence in regions where the influence stemss not only from spurious correlations but also from an actual physical influence of x; on Zz. Given a reasonable choice of system identification or model learning technique and a sufficiently rich excitation signal, this assumption will typically be satisfied in practice, as shown in the example below. It further follows that simulated trajectories in regions of local non-causality have a smaller MMD than trajectories in other regions.
Example (cont). Given assumption [ if we simulate experiments using equation 2o and compute the MMD for the resulting , the MMD will be lower if we choose ân (0) = 0 than for any other choice of ât (0). To volidate assumption 4 we identify the system using a regularization based nonlinear system identification algorithm called sporse identification of nontinear dynamics (SINDy) [Branton et al. (2016). By exciting the system for 500 time steps and using the SINDy implementation from de Silva et al. 2020)), we receive model. Simalating the MMD for different initial conditions (0) for (0) = -10 and x22 (0) = 10 then reveals that the MMD has a minimum at x (0) = 0 (cf- figure 3), i.e., the identified model satisfies assumption ] We will introduce the finite sample approzimation for the MMD used in figure D in section [4.4] and we provide a further ezomple empirically verify assumption 1 for a system with a hysteresis in appendiz D
We now specify the experimental design. To avoid regions of local non-causality, we propose to maximize the MMD given the model estimate. Thus, for checking whether æj causes 2i, choose input trajectories
In general, there may exist different. x for each combination x; and zy (and likewise for Une). We can also cover this case. However, we omit it here to simplify notation.
7

---

Published in Transations on Machine Learning Research (07 (07/2022)
-10-3
3
MMD
10
-8
-6
-2
2
x1
Figure 2: MMD of simulated experiments with different initial conditions z1) for the system from exam- ple 1 The MMD has a miriimmm at x1 (0) i.e., it reflects the local non-causalily. The fonile sample 2 approzimation MMD of the MMD used in this figure will be introduced in section 4.4
of length and initial conditions that solve
max (O), u2 an MMD (A. !!
(5a)
subject to (0) = (0) for all é f ul(() = (1) for all t,
We will discuss how to handle this optimization problem in practice section 5.2 In particular, in practice we typically do not require a global solution to equation 5a If we want check whether uy causes 2,, we choose input trajectories and initial conditions by solving
max (O) a0 II MMD(E., S! ()
(5b)
subject to (0) = z (0) for all w)(() = uf (t) for all + j,t.
Before stating our main theorem, we need further assumption. Given the system model from equation 1 could also have systems for which the influence of z, on z; becomes only apparent after a certain number
of time steps. To guarantee that identify the true causal structure, need to assume that this delay is at most equal to the length of an experiment T.
Assumption 2. Consider pair (z,,2;) for which æ, has cousal inffuence as per definition We assume that changes of æj, when outside of regions of local non-cauee non-causality, cause a change in after at most T discrete time steps.
We can now state the main theorem:
Theorem 1. Consider a completely e-controllable system described by equation with c = that fulfills assumptions - and 2 Let experiments be designed according to equation 3 for a fixed experiment length T and repeated infinstely often (í.e., we have a,, at! E RRxx with m >> oo). Then: MMD(a), a!1) if, and only if, zz respectively uj does not couse xi per definition 7
Proof. Let variables be non-causal. Then, have by definition 1 P(,; do(z)(0) =²(0))) P(::
do(x,(0)= = (0)). That is, the distribution of x, in both experiments is equal. Thus, MMD(x) = follows from Gretton et al. (012). Now, assume MMD(x,, = 0. This implics that the distribution of z; is equal in both experiments (Gretton et al. 2012), i.e., PP(x, | do(x,, (0) dt (0))) P(r; do(x,
(0)). This could be the case because (1) I3 and Ty are non-causal or (ii we have that x E Xnc and u E Unc. Due
to assumption [ and equation 5 there exists a for which either & Xne or u & Une Thus, we are in case (i) and r, does cause 24. The proof for uy and I; follows analogously.
Remark 3. Assumption 7 ensures that, if we design experiments using the optimézation algorithm equa- tion 5 we are not inside region of local non- causality. However, even if we weaken the assumption, the claim of theorem 7 would still hold. Since equation 5 finds the experiment design that mazimizes the MMD,
8

---

Published in Transations on Machine Learning Research (07/2022)
the algorithm fails if our model estimate assumes the influence of one variable on the other to be strongest in regions where there is actually no causal influence. Given reasomable model class and excilation signal for the initial model estimation, this should not occur in practice. We still keep assumption I in its stronger form as it allows us to also make causal inference in practical settings where we may not be able to globally solve equalion 5
Remark 4. Instead of doing mealtiple experiments of fixed length, one might consider doing one long exper- iment. This design has several disadvantages. First, if the process is stable, the probability distribution can become independent of the inilial conditions over time. Then, causal influences are only visible during the transient. Second, if we deal with non-ergodic systems, where time and spatial overage are not the some, we only have a valid testing procedure if we do maltiple runs of each experiment. Lastly, doing multiple experiments ensures that the different runs are independent of each other. Data generated by single long experiment can be correlated, and oblaining a valid test is way more involved (Solowjow el al. 2020).
4.2 e-Controllability
For a stochastic system as in equation [ it is in general impossible to steer the system exactly to the initial conditions suggested by equation 5 i.c., we need to resort to controllability with > O (cf. dcfinition 2 Nevertheless, even in such cases, it is still possible to guaranlee the consistency of the causality testing procedure, However, we need a stricter definition of controllability
Definition 4. Let the system given in equation 1 be e-controllable according to definition g and consider some arbitrary, but fixed di, Then, the system given in equation 7 is said to be completely e-controllable in distribution if, for any z²) and any Zdes with De,des = =deee there erists an imput sequence u(0:t) such that se((e) always follows the same distribution; i.e., (P( = TP* for some TP" that does not depend on 200) or any component of Ides except Lt,des = =iides and Pr ([(]) - Ides l12 2 c) 1 111 with 0 < < 1.
In other words, the definition states that, for any z0)), we can generate input trajectories that guarantee that the fixed component Fi of Fdee is matched in distribution. Linear systems with additive Gaussian noise that are controllable following Kalman (960a)) are also controllable in the sense of definition 4 as we show in the appendix. We further need to make an assumption about the initial conditions suggested by equation 5 To steer the system to those initial conditions without ending up in a region of local non-causality, we need that states that are c-close to those initial conditions are not in a region of local non-causality.
Assumption 3. For the initial conditions suggested by equation 5 we assume that there exist open balls with radius Ve centered around each element of (0) and (0) that are outside of the local nor- non-causality sets.
This assumption is not very strong since equation 5 suggests initial conditions for which the influence of Ty respectively uy on x, is maximal. Thus, it is unlikely that these initial conditions are close regions of local non-causality. The main reason why this assumption is needed is to exclude corner cases in which the causal influence only exists in single points of the state and input space. In such cases, the probability of successfully stecring the system to those points is zero. However, considering the variables as non-causal may then anyway be reasonable. Equipped with these three assumptions, we can now state:
Corollary 1. Consider system as in equation 7 that is completely e-controllable in distribution according to definition D and fulfills assumptions 1 g and 3 Let experiments be designed as in equation 5 for faed experiment length T, trajectories that steer the system to the initial conditions of the experiments be chosen such that (0)) = P(x (0) for all E + j. and experiments be repeated infimitely often. Then: MMD(a), a!!) 0 if, and only if, 2j respectively uj does not cause I, according to Definition
Proof. Let variables be non-causal. Then, we have P(x((0)) = P(x! for all E + j. thus, also the distribution of the obtained x; trajectories is equal and we have MMD(x(, = 0 (Gretton el al. 2012). Now, assume xIM((x..!") = O. This implies that the distribution of x; is equal in both experiments (Gretton et al. 2012). By assumption | existing local non-causalitiess are reflected by the model and thus, equation 5 will suggest experiments outside of such regions. Assumption ensures that we can steer the system to those regions. Thus, if distributions are equal, non-causality must be global in definition 1 7
9

---

Published in Transations on Machine Learning Research (07/2022)
4.3 Linear Systems
Local non-causality in definition is nonlinear phenomenon. If we assume equation 1 lo be linear time- invariant (LTI) with Gaussian noise u(t), we can reveal the true causal structure without the optimization procedure in equation 5 making this case substantially easier. For an LTT system, equation | reads
x(t) - A'x(0) Σ A -(t---1
(6)
with state transition matrix A R input matrix B E Raxm and v(t) (O Σ.). The system in equation [6 is controllable per definition 4f it satisfies the classical controllabiltyy condition from Kalman (1960a), Le., if the matrix (B AB A)-1B) has full row rank, as we show in lemma in the appendix. We can then state the following theorem, whose proof is provided in the appendix:
Theorem 2. Assume an LTI system as in equation G whose A and B matrices satisfy Kalman's con- trollability condition. Let experiments be designed as in equation 2a and equation 2b respectively. Then: MMD(2 a = 0 if, and only if, zj respectively uj does not couse Vs as per definition
4.4
Test with Finite Samples
Until here, we derived guarantees in the infinite sample limit. In practice, we can only carry out finitely many experiments, ie., we have finitely many samples of the random variable sequence z, Thus, we need a finite sample approximation of the MMD.
Lemma 1. Consider m experiments with fixed length T, i.e., a E RRxT but now we also have m < co. An unbiased empirical estimate of the squared population MMD can be computed as
MMD (a) x m(m 1)
Σ (k(" " k("a a $ *x A(T *z!1) - k("at, xz))))
(7)
where "a! denotes element of the x, trajectories from experiment I.
Proof. The m trajectories for xt and all are independent and follow pl and pd respectively, Thus, we are in the same setting in (Gretton et. al. 2012 lem. 6) and the proof follows shown therein.
7
For a finite sample approximation, we, in general, have MMD (x}, x!') > even if plp-= pll Thus, we need do a hypothesis test and derive test statistic. We assume the mull hypothesis
Ho: p'=
(8)
and obtain the test statistic from the following result. Lemma 2. Assume 0 < k("al, *xfl) < K. Then
->2m2
rr,f ail
( (MMD* (2), a MMD? " ,Ill) > 7 S
exp
8K2
where ma = [m/2]. The hypothesis test of level for the null hypothesis in equation 8 has the acceptance region
(+m Vlog(a-1))
4K
MMD (2), a!") <
Proof. The trajectories for x and T are independent and follow pl and p" respectively. Thus, the proof follows from theorem 10 and corollary 11 in Gretton et al. (2012).
Note that the unbiased estimate lemma can even become negative, the discussion after lemma Gretton
2012
10

---

Published in Transations on Machine Learning Research (07/202)
5 Implementation
The results in section show that we are able to detect whether the variables Ij or Wj have a causal influence on Ty. In practical implementations, two challenges remain: first, we want to minimize the number of experiments we need to carry out on the physical platform. Second, we may be unable to obtain a global solution to equation
5.1 Heuristic Test with Finite Samples
The test statistic provided in lemma 2 enjooy a theoretical foundation. However, the threshold decreases as m -///,, i.c., we need many exporiments to not be overly conservative. While more efficient test statistics exist (e.g,; achieved through subsampling discussed in section 6 of Gretton al. (2012)), generating all this data through experiments on real-world systems is often undesirable, e.g., because it is time-consuming and may cause excessive wear and tear on the hardware. Thus, propose an alternative test statistic that can be obtained from the model estimate f. This alternative test statistic is heuristic and forgoes the theoretical properties but is efficient to implement and yields good results in practice, as we show in section G
We estimate a model fine that assumess x; and 2j respectively "j are non-causal (i.e., we do not use the data Ij respectively uj when estimating fr.n)) We propose to use this model to decide whether to accept the null hypothesis of 24 and 2j respectively u, being non-causal. That is, we replace our current model fo for state component oj with fine
2 MMD a <E|MD af! fe,no)] + vy Var MMMD (e), 2ll fi,no))-
2
(9)
Expected value and variance in equation 9 can be estimated through Monte Carlo simulations. For these simulations, we use the true initial conditions x¹ (0) and z²))). That way, we account for uncertainty due unequal initial conditions between experiments. The signifieance level of the test be adjusted through using Chebyshev's incquality (Chebyshev 1867
5.2 Experimental Design
The framework for testing causality between state variables is summarized in algorithm 1 (and works analo- gously for inputs). After having obtained an initial model (see 11. 1-2 in algorithm we run equation 5a to find initial conditions and input trajectories for testing non-causality of one specific combination of ze and zj (I. 6). The optimization in equation 5a may be arbitrarily complex or even intractable, depending on the chosen model class. However, finding a global optimum of equation 5a is not necessary. The goal of the opti- mization procedure is to avoid regions of local non-causality. We thus optimize equation 5a until is above a threshold & to be confident that we are not in a region of local non-causality. In practical applications, we can often already achieve this through a reasonable initialization of the optimization problem, i.e., by choosing initial conditions for Hy far apart possible. We then run the designed experiment and collect the data (1. 7). Tdeally, we would like to use data from this single experiment to test for the causal influence of z on all other state components. Thus, we check for which the experiment yiclds an expected MMD that is higher than a second threshold 62 and do the hypothesis lest for all of those (II. 9-13).
6 Evaluation
We evaluate the framework on three systems. First, we identify the causal structure of one arm of the robot Apollo (Kppleer et al.) (018)) shown in figure 3 in section 6.1 Then, we demonstrate the causal identification of simulated quadruple tank process (ef. figure 5) in section 6.2 In both cases, we present the general setup and results in the following, while we defer implementation details to the appendíx. Lastly, we discuss synthetic linear toy example. We use this example to highlight once more the importance of considering notion of controllability and to compare our method to sparse identification and a causal discovery algorithm For all experiments, we use a Gaussian kernel to compute the MMD.
3Code for both simulation examples available http:r///fttii.ccaaaammmmmaaaait causal hyt.eryyyyyyyaaaaaaeettttturrrrr sttuccure
11

---

Published in Transations on Machine Learning Research (07/2022)
Algorithm Pseudocode of the proposed framework.
Excite system with imput signal, collect data
2: Obtain f through black-box system identification for Iy in do
4: -Itest = [fI,...,
> states to be tested for non-causality
for my Steet do
design experiment test whether æj causes
6:
Run equation 5a until E[MMD (#1, [f) > 61 Run causal experiments, collect data
8:
for 24 in Ptest do
The experimental designed for 2, and ze might also yield a valid test other i > If the empirical MMD is above the threshold for some other 2i, also test that 24
9:
TEMMM (#), a! f) > 82 then
10:
Obtain fine
11:
Obtain test statistic via Monte Carlo simulations if equation
12:
g holds then
independence test
13:
f. frac
14:
Delete x; from Zeest
Figure 3: The robot showing initial postures for two experiments.
6.1 Robot Experiments
We consider kinematic control of the robot; that is, we can command desired angular velocitics to the joints, which are then tracked by low-level controllers (taking care, among other things, of the robot dynamics) [Berenz et al. 2021). As measurements, we receive the joint angles. The goal of the causal identification is to learn which joints influence each other and which joints inffuenced by which inputs. We consider four joints of the robot arm in the following experiments. From the robot arm's design, we know its kinematic structure, which is described by Φ; M4 for each joint, where Φi denotes the angle of joint i and the control input. That is, we expect each joint velocity di to depend only on the local input u, and not other variables. While the dynamics are approximately linear, we do not rely on this information and are thus in the setting discussed in section 4.2 the following, we will investigate whether our proposed causal identification can automatically reveal this structure.
Following algorithm 1 we start by identifying model f. In this experiment, estimate linear state- space model. As expected, the initial model suggests that all joints are linked to each other and all inputs due to spurious corrolations. We then design experiments for causality testing, example trajectorics of such experiments are shown in figure (left). The empirical squared MMD of the resulting trajectories is compared with the test statistic. The trajectories in figure | (left) already suggest that the experiments are in line with the kinematic model: while the two trajectories of joint for different initial conditions of joint are essentially equal (blue dashed and green dotted lines overlap), the trajectories of joint for different choices of the third input fundamentally different. This is also revealed through the proposed causality
12

---

Published in Transations on Machine Learning Research (07/2022)
test. The middle plots of figure show the empirical squared MMD (left-hand side of equation 9 and the test threshold (right-hand side of equation [] for the experiments that were conducted to test the influence of the initial conditions of the third joint (top) and of the third input signal (bottom) on all joints. As can be seen, the causal identification reveals that the third joint does not inffuence any other joint, and the third input only affects the third joint. Note that the third joint's trajectories are obviously different when choosing different initial conditions for the third joint. However, since this expected, we subtract the initial condition in this case to investigate whether the movement starting from these distinet initial conditions differs. The remaining experiments (results are contained in the appendix) yield similar results. In summary, the causal identification successfully reveals the expected causal structure.
The experiment design in equation 5 requires us to solve an optimization problem. Nevertheless, algorithm [ introduces two parameters, 81 and 52. These can be used to stop the optimization early in case the predicted MMD high enough for us to be confident that we are not a region of local non-causality (6i) and that we can use the design to test for all influences of joint or input signal (6z). To design 81 and 82, we can look at the system's 's noise level and choose them some orders of magnitude higher. As discussed in section [ high predicted MMD can often already be achieved through sophisticated guesses, i.e., by choosing initial conditions far apart from each other and diverse input signals. We follow this approach. When testing for the influence of the third joint on others and choosing initial conditions that are far apart, we predict MMDs of around 0.5. The model we estimate initially for the robot arm has a noise standard deviation below 1 10-4, That is, the predicted MMD is way above the noise level of the system and way above the MMDs we find in experiments. Thus, for any choice of d1 and 62 below 0.5, we can confidently accept this experiment design. If we were even more conservative and chose them above 0.5, we would need to optimize the experiment design further.
To investigate the generalization capability, we compare predictions of the model finit obtained from the initial system identification and the model feaus that was learned after revealing the causal structure. We use the same training data to estimate the model parameters in both cases. However, for feany, we leverage the obtained knowledge of the causal structure when estimating parameters. In contrast, for finit. we do not take any prior knowledge into account. As test data, we use an experiment that was conducted to investigate the influence of the initial condition of joint 3 the other joints and let both models predict the trajectory of joint 1. For this experiment, the initial angle of joint 3 is close to its maximum value, case that is not contained in the training data. As can be seen in figure 4 (right), the predictions of feaus (blue) are very close to the true data (green, dashed), i.c., the model can generalize well, while the predictions of fini deviate significantly.
6.2 Quadruple Tank Process
The experimental demonstration in the previous section showed that the presented algorithm can successfully identify a real-world robotic system's causal structure. However, the causal structure is relatively straight- forward, and dynamics are approximately linear. To stress the method's ability to perform similarly well more complex structures and with nonlinear dynamics, we now consider the quadruple tank system from Johansson 2000), illustrated in figure 5 Its contimuous-time dynamics are given by
41
V /2gm
a3
V /2gugg +
A1 W1
S1
A
A
02
a3
Ç2 U2
22 =
V 2gxz
+
/2g01
+
A
Ag
A2
(10)
43
1 Co
43
v 22gz
Ag
Ag
04
1 1-6 C1
2gr
U1.
A
A₁
where x; denotes the water level of tank i, u, the flow rate of pump i, g the gravitational constant (9.81 ms= The dynamics of the quadruple tank process are, thus, clearly nonlinear. That is, we cannot expect good performance if we approximate them using linear state-space
and the remaining constants as in figure 5
13

---

Published in Transations on Machine Learning Research (07/2022
Trajectories from causality experiment
Causality test
Predictions with learned model
Angie Joint (rad)
1.65
5.3e-T
Ángie Joint I (rad)
MMD
1.6
1e-16
20
1.00
-1e-16
40
20
40
60
80
20
10
60
80
100
Joint
Angle Joiat 1 (rad)
109
Augle Joint 1 [rad]
1.7
1.65
Exp.
MMD
Exp. II
Exp. MMD
1.6
Test stat.
Exp.
-0.1
1.55
Exp. 11
10-14
20
40
60
80
20
40
60
80
100
# Joint
Figure Causality tests and model evaluation for the robotic system. Plots on the left show example trajectories of two experiments, in the middle the experimental MMD and the test threshold for joint and on the vight predictions based on the iniläal model and on the refined model after the causal idenlifieation.
A3
A
Tank 3
Tank 4
C1
A
Aa
Tank
Tank 2
Figure 5: Schematic of the quadruple tank process from Johansson (2000).
model. Therefore, discretize the system and use Gaussian processes (GPs) to model the dynamics, Williams & Rasmussen (2006) for an introduction and details. In particular, we model each state of the system with a GP, where each GP uses all states and inputs to predict its assigned state variable at the next time step.
We procoed similarly as for the robot experiments. We select initial conditions through sophisticated guesses, i.c., far apart from each other. In simulation, we can initial conditions and do not need to steer the system there. However, be more realistic, sample initial conditions from normal distribution
14

---

Published in Transations on Machine Learning Research (07/2022)
10a
10
Exp. MMD
MMD
Test stat.
MMD
10
10-4
10
10
2
4
2
3
# Tank
# Tank
Figure 6: Causality tests for the quadruple tank system. The plots show the results of the experiments testing the influence of the third tank on all others (left) and the influence of un on all tanks (right).
with mean equal to the selected initial conditions and variance of 10-2. Thus, we only consider c- controllability.
Results of the testing procedure are shown in figure G The left plot shows test statistic and experimental MMD for the third tank's influence on all others. It reveals that the third tank influences itself and the first tank, which is in line with the schematic in figure 5 and the dynamics in equation [0] The right plot illustrates the test statistic and experimental MMD for the influence of 41 on the tanks. Also these results are accordance with figure 5 as the experiments suggest that u1 influences all but the third tank. The remaining results are collected section [ in the appendix. Similarly as the ones in figure they reveal the causal structure that can be inferred from figure
In this case, the causal identification results in reduction of computational complexity. Especially for the third and the fourth tanks, which are only influenced by U2 respectively and by themselves, the input dimension of their corresponding GPs decreases from 10 to 3. The complexity of standard GP regression grows with Onx²) with the number of datapoints n. Thus, if we can reduce the dimensions that need to be considcred and, with that, the number of data points that we use for the regression by 70 % as in this example, we considerably reduce the computational complexity.
6.3 Synthetic Example and Comparison
Lastly, we present a synthetic, linear example. We consider an LTI system in equation 6 with
0.9
-0.75
1.2
0.00
0
0
A = 0
0.9
-1.1
B =
0
0.06
0
(11)
0
0.7
0.07
0.05
and Gaussian noise with standard deviation 1 x 10-4, For this example, we apply algorithm without the need for the optimization procedure since the example is linear. Again, we want to stress the importance of an appropriate notion of controllability. That is, instead of assuming that we can set the system to initial conditions, we always steer the system the initial conditions required for that experiment. For this, we employ an approach set-point tracking that has, for instance, been discussed in Pannocchia et al.) (2005). Given a desired state Ides; seek feedback control law of the form that depends both the desired state and the current state. We obtain the gain matrix F' using standard methods from linear optimal control (Andrson & Moore 2007), in particular, the linear quadratic regulator (LQR). Thus, we can rewrite the incremental dynamics of the system as
Mr + Fxr, i.e., control law
x(1 11 = Ax(() + BM + v()),
(12)
where À == A + BF. We choose the feodforward term M such that the reference is matched in stationarity, i.c., want to achieve
x (1 - Ba
(13)
15

---

Published in Transations on Machine Learning Research (07 (07/2022)
2
Causal SINDy
0
20 40 60
80 100
Figure 7: Comparison prediction capabilities. In blue predictions of the true model obtained after causal identification, in green the model obtained from SINDy with parameter choices from Brunton et al. (2016). We do not compare to PCMCI in the numerical experiment since PCMCI only reveals which causal influences exist, but not how strong they ere, i.e., different from algorithm and SINDy, it does not return full system model that can be used in simalations.
Thus, we need
M ((I A)-
(14)
to track the reference point. To compute M, we use the matrices À and B of the estimated model f. We start the experiment once ||x(t) - -deslll < 0.01.
Analyzing the resulting data with the MMD lets us, similar as before, infer the true causal structure. particular, the causal analysis reveals that æj does not cause z2 nor I3, 22 does not cause T3, and U2 does not cause 23.
In addition to those results, compare our method with a sparse identification method and causal discovery algorithm for this linear example. In both cases, we excite the system for 100 000 time-steps with chirp as an input signal.
First, we use SINDy, as in example to identify the underlying discrete-time system Sparse identification algorithms seek to find a trade-off between model complexity and goodness of fit. Thus, for such algorithms, parameter needs to be selected to indicate how important accuracy is compared to complexity. When applying SINDy the synthetic system, its success in finding the true causal structure depends on the choice of this parameter. We start with the parameter scttings that were used for thrce-dimensional lincar example in Brunton et al. (2016). For those setlings, the algorithm does not recognize the influence of u on T1, ie., the first row of the B matrix consists of zeros. In figure 7 we show the effects that this error can have by comparing predictions of the true model obtained after the causal identification with predictions of the model obtained from SINDy. For this, we set 41 = 10, U2 = 0, and Ug = - 14 and simulate both models for 100 time steps. As the SINDy model does not reflect the influence of リ₁ on 31, assumes that Ii does not move. However, the model obtained after the causal identification reflects this influence and, thus, correctly predicts the movement of 21.
Only when lowering the threshold parameter can SINDy recover the true causal structure of the system. This stresses the general shortcoming of sparse identification methods when identifying the causal structure of control system. Depending on the parameter set settings, they may or may not recover the true causal structure. However, their general purpose different. They seek a trade-off between model complexity and accuracy. Thus, neglecting a causal link with comparably weak inffuence might be the desired outcome. In contrast, we seek the true causal structure, independent of how strong the influence is.
We the implementation provided in de Silva
2020
16

---

Published in Transations on Machine Learning Research (07/202)
Second, we compare our algorithm to the PCMC algorithm proposed in Runge et al 2019b 6 algorithm focuses on detecting causal influences from time-series data. That is, it does not yield an identified system matrix but discovers which variables have a causal influence on which others. When again exciting the linear system with a chirp signal and running the PCMCI algorithm, we receive an output that we can interpret as
This
11 (k- + 1) = a1 (k) + (2) + aza3 (k) + b1 (1
22 + 1) = 01 (k: 1) + 05 (k) + 25 (k) + ak 22 2) + b2u2(k) + byus(k -1 1)
(15)
ra(k + 1) = u7 (A: 4) + agrz(k) (k) + ( + (k))
where all weights a, and 6; are non-zero. Besides, the algorithm discovers an influence of on Ug. However, since such links were ruled out by design for the other algorithms, we neglect this here. Nevertheless, also the equations as stated above are not in line with the actual matrices in equation IT For instance, following equation 15 the variable 21 has a causal influence on 22 and 23 through the non-zero factors 04 and a7. However, given the upper-triangular structure of the A matrix, 23 not influenced by any other state variable apart from itself. Similarly, also d4 should be zero. Thus, the algorithm does not reveal the true structure of the system.
7 Conclusion
We presented method that identifies the causal structure of dynamical control systems by conducting experiments and analyzing the generated data with MMD-based techniques. It differs from prior approaches to causal inference in that it uses a controllability notion that is suitable to design experiments for control systems. We evaluated the method on a real-world robotic system and simulated quadruple tank system. Our algorithm successfully identified the underlying causal structure of both systems, which in turn allowod to learn model that accurately generalizes to previously unseen states while redueing computational complexity.
Acknowledgments
The authors would like to thank Manuel Wüthrich and Alonso Marco Valle for insightful discussions, Steve Heim for valuable feedback, and Vincent Berenz for support with the robot experiments. This work has received funding from the German Research Foundation within the SPP 1914 (grant TR 1433/1-1), the Federal Ministry of Education and Rescarch (BMBF) and the Ministry of Culture and Science of the German State of North Rhine-Wierr Westphalia (MKW) under the Excellence Strategy of the Federal Government and the Länder, the Cyber Valley Initiative, the Max Planck Society, the Knut and Alice Wallenberg Foundation, and the Swedish Research Council.
References
Hirotogu Akaike. Information theory and an extension of the maximum likelihood principle. In International Symposium on Information Theory, pp. 267-281, 1973.
Brian DO Anderson and John B Moore. Optimal Control: Linear Quadratic Methods. Courier Corporation, 2007.
Agamirza E Bashirov and Kerim R Kerimov. On controllability conception for stochastic systems. SIAM Journal on Control and Oplimization, 35(2):384 398, 1997.
Vincent Beronz, Maximilien Navcau, Felix Widmaier, Mamuel Wüthrich, Jean-Claude Passy, Simon Guist, and Dieter Büchler. The o80 C++ templated toolbox: Designing customized Python APIs for synchro- nizing realtime processes. Journal of Open Source Software, 6(66):2752, 2021.
an adaptation of the PC algorithm, uamed after its creators Peter Spirtes and Clark Glymour Spirtes al. 2000 which momentary conditional independence (MCI) test deal with auto-correlations in time-scries data For our simulations, use the implementation and parameter settings provided https: http://ggttttbttt/ com/ jakobrunge/ tigranite
17

---

Published in Transations on Machine Learning Research (07/202)
Stefano Boccaletti, Vito Latora, Yamir Moreno, Martin Chavez, and Dong-Uk Bnng-Uk/ Hwang. Complex networks: Structure and dynamics. Physics Reports, 424(4-5):175-38) 308, 2006.
Josh Bongard and Hod Lipson. Automated reverse engineering of nonlinear dynamical systems. Proceedings of the National Academy of Sciences, 104(24):9943 9948, 2007.
Steven L Brunton, Joshua Proctor, and J Nathan Kutz. Discovering governing equations from data by sparse identification of nonlinear dynamical systems. Proceedings of the National Academy of Sciences, 113(15):3932 3937, 2016.
Pafiutti Lvovich Chebyshev. Des valcurs moyennes. Journal de Mathématiques Pures et Appliquées, 2(12): 177 184, 1867.
Zhitang Chen, Kun Zhang, Laiwan Chan, and Bernhard Schölkopf. Causal discovery via reproducing kernel Hilbert space embeddings. Newral Compulation, 26(7):1484-1511 2014.
Kacper Chwialkowski and Arthur Gretton. A kernel independence test for random processes. In International Conference Machine Learning, pp. 1422- 122--430, 2014.
Kacper Chwialkowski, Dino Sejdinovic, and Arthur Gretton. A wild bootstrap for degenerate kernel tests. 2014.
Tn Advances in Newral Information. Processing Systems, pp.
Brian de Silva, Kathleen Champion, Markus Quade, Jean-Christophe Loiseau, J. Kutz, and Steven Brunton. PySINDy: A Python package for the sparse identification of nonlinear dynamical systems from data. Journal of Open Source Software, 5(49):2104, 2020.
Selva Demiralp and Kevin D Hoover. Searching for the causal structure of a vector autoregression. Oxford Bulletin of Economics and Statistics, 65:745 767, 2003.
Michael Eichler. Graphical Gaussian modelling of multivariate time series with latent variables. In Interna- tional Conference on Artificial Intelligence and Statistics, pp. 193-200 2010.
Michael Eichler. Causal inference in time series analysis. In Carlo Berzuini, Philip Dawid, and Luisa Bernardinell (eds.), Causality: Statistical Perspectives and Applications, chapter 22, pp. 327-354. John Wiley & Sons, Ltd, 2012.
Abbas Emami-Naeini and G Franklin. Deadbeat control and tracking of discrete-time systems. IEEE Transactions on Automatic Control, 27 1982.
Doris Entner and Patrik O Hoyer. On causal discovery from time series data using FCL. Probabilistic Graphical Models, pp. 121-128, 2010.
Karl J Friston, Lee Harrison, and Will Penny. Dynamic causal modelling. Neuroimage, 19(4):1273 1302, 2003.
Kenji Fukumizu, Arthur Gretton, Xiaohai Sun, and Bernhard Schölkopf. Kernel measures of conditional dependence. In Advances in Neural Information Processing Systems, pp. 489-496, 2008.
Arthur Gretton, Karsten M Borgwardt, Malte J Rasch, Bernhard Schölkopf, and Alexander Smola. kernel two-sample test. Journal of Machine Learning Research, 13(Mar):723 773, 2012.
Joao Hespanha. Linear Systems Theory. Princeton University Press, 2018.
Karl Henrik Johansson. The quadruple-tank process: A multivariable laboratory process with an adjustable zero. IEEE Transactions on Control Systems Technology, 8(3):456 06016646556666 2000.
Rudolf Emil Kalman. Contributions the theory of optimal control. Boletin de la Sociedad Matematica Mezicuna, 5(2):102-11) 119, 1960a.
Rudolf Emil Kalman. On the general theory of control systems. In International Conference on Automatic Control, pp. 1960b.
18

---

Published in Transations on Machine Learning Research (07/2022)
Daniel Kappler, Franziska Meier, Jan Issac, Jim Mainprice, Cristina Garcia Cifuentes, Manuel Wüthrich, Vincent Berenz, Stefan Schaal, Nathan Ratliff, and Jeannette Bohg. Real-time perception meets reactive motion generation. IEEE Robotics and Automation Letters, 3(3):1864 1871, 2018.
Solomon Kullback and Richard A Leibler. On information and sufficiency. Statistics, 22(1):79 86, 1951.
The Annals of Mathemalical
Fabien Lauer and Gérard Bloch. Hybrid System Identification: Theory and Algorithms for Learning Switching Models. Springer, 2018
Hui Liu, Jun-An Lu, Jinhu Lü, and David J Hill. Structure identification of uncertain general complex dynamical networks with time delay. Automatica, 45(8):1799 1807, 2009.
Lennart Ljung. System Identification: Theory for the User. Prentice Hall PTR, 1999.
David Lopez-Paz, Krikamol Muandet, Bernhard Schölkopf, and llya Tolstikhin. Towards learning theory of cause-effect inference. In International Conference Machine Learning, pp. 1452-1461, 2015.
Daniel Malinsky and Peter Spirtes. Causal structure learning from multivariate time series settings with unmeasured confounding. In ACM SIGKDD Workshop on Causal Discovery, pp. 23-47, 2018.
Donatello Materassi and Giacomo Innocenti. Topological identification in networks of dynamical systems. IEEE Transactions on Automatic Control, 55(8):1860-188) 1871, 2010.
Alessio Moneta, Nadine ChlaB, Doris Entner, and Patrik Hoyer. Causal search in structural vector autore- gressive models. In NIPS Mini-Symposium on Causality in Time Series, pp. 95-1944 2011.
Joris M. Mooij, Dominik Janzing, and Bernhard Schölkopf. From ordinary differential equations to structural causal models: The deterministic case. In Conference on Uncertainty in Artificial Intelligence, pp. 440 040-448 2013.
Duy Sguyen-Tuone and Jan Peters. Model learning for robot control: A survey. Cognétive Processing, 12 (4):319-340 2011.
Gabriele Pannocchia, Nabil Laachi, and James B Rawlings. candidate to replace PID control: SISO- constrained LQ control. AICHE Journal, 51(4):1178-11) 11888,111,,1110 2005.
Judea Pearl. Causal diagrams for empirical research. Biometrika, 82(4):669 688, 1995.
Judea Pearl. Theorctical impediments to machine learning with seven sparks from the causal revolution. arxiv preprint ar arXic:1801.00016 2018.
Judea Pearl and Dana Mackonzie. The Book of Why: The New Science of Cause and Effect. Basic Books, 2018.
Jonas Peters, Dominik Janzing, and Bernhard Schölkopf. Causal inforence on time scries using restricted structural equation models. Advances in Neural Information Processing Systems, pp. 154 162, 2013.
Jonas Peters, Dominik Janzing, and Bernhard Schölkopf. Elements of Causal Inference: Foundations and Learning Algorithms. MIT Press, 2017.
Jonas Peters, Stefan Bauer, and Niklas Pfister. Causal models for dynamical systems. In Hector Geffner, Rina Dechter, and Joseph Y. Halpern (eds.), Probabilistic and Causal Inference: The Works of Judea Pearl, pp. 671-690. 2022.
Niklas Pfister, Stefan Bauer, and Jonas Peters. Learning stable and predictive structures in kinetic systems. Proceedings of the National Academy of Sciences, 116(51))25405-25411 2019.
Christopher J Quinn, Todd P Coleman, Negar Kiyavash, and Nicholas G Hatsopoulos. Estimating the directed information to infer causal relationships in ensemble neural spike train recordings. Journal of Compulational Neuroscience, 30(1):17-44 2011.
19

---

Published in Transations on Machine Learning Research (07/2022)
Jacob Roll, Alberto Bemporad, and Lennart Ljung. Identification of piecewise affine systems via mixed- integer programming. Aulomalica, 40(1):37-50; 2004.
Paul K. Rubenstein, Stephan Bongers, Bernhard Schölkopf, and Joris M. Mooij. From deterministic ODEs dynamic structural causal models. In Conference on Uncertainly in Artificial Intelligence, 2018.
Cynthia Rudin. Stop explaining black box machine learning models for high stakes decisions and use inter- pretable models instead. Nature Machine Intelligence, 155):206 215, 2019.
Jakob Runge, Sebastian Bathiany, Erik Bollt, Gustan Camps-Vaie Dim Coumou, Ethan Deyle, Clark Gly- mour, Marlene Kretschmer, Miguel D Mahecha, Jordi Muñoz-Mari, Egbert Nes, Jonas Peters, Rick Quax, Markus Reichstein, Marten Scheffer, Bernhard Schölkopf, Peter Spirtes, George Sugihara, Jie Sun, and Jako Zscheischler. Inferring causation from time series in Earth system sciences. Nature Communications, 10(1):1-13, 2019a.
Jakob Runge, Peer Nowack, Marlene Kretschmer, Seth Flaxman, and Dino Sejdinovic. Detecting and quanti- fying causal associations in large nonlinear time series datasets. Science Advonces, 5(1 ((11))eaa44996, 2019b.
Cristopher Salvi, Maud Lemercier, Chong Liu, Blanka Horvath, Theodoros Damoulas, and Terry Lyons. Higher order kernel mcan embeddings to capture filtrations of stochastic processes. Advances in Neural Information Processing Systems, 2021.
Michael Schmidt and Hod Lipson. Distilling free-form natural laws from experimental data. Science, 324 (5923):81-88) 85, 2009.
Bernhard Schölkopf. Causality for machine learning. In Hector Geffner, Rina. Dechter, and Joseph Y. Halpern (ods.), Probabilistic and Causal Inference: The Works of Judea Pearl, pp. 765-804. 2022.
Johan Schoukens and Lemnart Ljung. Nonlinear system identification: A user-oriented road map. 1EEE Control Systems Magazine, 39(6)228-99, 2019.
Gideon Schwarz. Estimating the dimension of model. The Annals of Statistics, 6(2):461-466 1978.
Shahin Shahrampour and Victor M Preciado. Topology identification of directed dynamical networks via power spectral analysis. IEEE Transactions on Automatic Control, 60(8):2260-222) 2265,,2225522220 2014.
Karthikeyan Shanmugam, Murat Kocaoglu, Alexandros G Dimakis, and Sriram Vishwanath. Learning causal graphs with small interventions. Advances in Neural Information Processing Systems, pp. 3195-3203,
2015.
Max Simchowitz, Horia Mania, Stephen Tu, Michael I. Jordan, and Benjamin Recht. Learning without mixing: Towards a sharp analysis of linear system identification. In Conference on Learning Theory, pp. 439-473, 2018.
Alexander Sokol and Niels Richard Hansen. Causal interpretation of stochastic differential equations. Elec- tronic Journal of Probability, 19(100):1-24, 2014.
Friodrich Solowjow, Dominik Baumann, Christian Fiedler, Andreas Jocham, Thomas Scel, and Sebastian Trimpe. A kernel (wo-sample test. for dynamical systems. ai preprint arXiv:2004.11 11098, 2020.
Peter Spirtes, Clark N Glymour, Richard Scheines, David Heckerman, Christopher Meek, Gregory Cooper, and Thomas Richardson. Cousation, Prediction, and Search. MIT Press, 2000.
Bharath K Sriperumbudur, Arthur Gretton, Kenji Fukumizu, Bernhard Schölkopf, and Gert RG Lanckriet. Hilbert space embeddings and metrics on probability measures. Journal of Machine Learning Research, 11(Apr):151-11661, 2010.
Bharath K Sriperumbudur, Kenji Fukumizu, and Gert RG Lanckrict. Universality, characteristic kernels and RKHS embedding of measures. Journal of Machine Learning Research, 12(7), 2011.
20

---

Published in Transations on Machine Learning Research (07/2022
Bharath K Sriperumbudur, Kenji Fukumizu, Arthur Gretton, Bernhard Schölkopf, and Gert RG Lanckriet. On the empirical estimation of integral probability metries. Electronic Journal of Slatistics, 6:1550-1599, 2012.
Klaas Enno Stephan, Will D Penny, Rosalyn Moran, Hanneke EM den Ouden, Jean Daunizeau, and Karl J Friston. Ten simple rules for dynamic causal modeling. Neuroimage, 49 40(4):3099-3109, 2010.
Yoshifumi Sunahara, Teruo Kabeuchi, Yoshiharu Asada, Shin Ichi Aihara, and Kiyotaka Kishino. On stochastic controllability for nonlinear systems. 1EEE Transactions Automatic Control, 19(1):49 54, 1974.
Zoltán Szabó and Bharath Sripcrumbudur. Characteristie and universal tensor product kernels. Journal of Machine Learning Research, 18:233, 2018.
Paul MJ den Hof, Arne Dankers, Peter SC Heuberger, and Xavier Bombois. Identification of dynamic models in complex networks with prediction error methods- basic methods for consistent module esti- mates. Automatica, 49(10))29943300) 2013.
Christopher KI Williams and Carl Edward Rasmussen, Gaussion Processes for Machine Learning, volume 2. MIT press Cambridge, MA, 2006,
Karren Yang, Abigail KatcofT, and Caroline Uhler. Characterizing and learning equivalence classes of causal DAGs under interventions. In International Conference Machine Learning, pp. 15411,550, 2018.
Dongehuan Yu. Estimating the topology of complex dynamical networks by steady state control: Generality and limitation. Automatica, 46(12):2035-2044 2040, 2010.
A Proof of Theorem 2
An LTI system with Gaussian noise follows a normal distribution, whose mean and variance are given by
B[x(()) = A'x(0) +
Σ
A' Bu(t 1-
(16a)
i=u
Var(x Σ
V'raffu(t - i)]
(16b)
where we assume Var[u(t)] = E, for all t. For such systems, we will first show that if equation 16a obeys the controllability conditions stated by Kalman, the system is also controllable according to definitions and
Lemma 3. The system in equation G is completely 6-controllable in distribution if the deterministic part obeys the controllability condition stated in Kalman (1960a).
Proof. The expected value equation 16a represents the deterministic part of the system. Thus, according to Kalman 1960a), we can design an input trajectory that steers equation 16a to any point in the state space. Since do not assume constraints on the input the state variables, the desired state can be reached with trajectory within do time steps (cf. deadbeat control (Kalman 1960b Emami--aeee &
Franklin 1982 That is, starting from any x00) E X' we can steer the system state toward x(g) _ x in time steps and obtain the distribution
E [x(()]] I
g-1 Var[z(g) = Σ A'E..
T=0
The probability ||=(g) 2²1l² being larger than e is given by the cumulative distribution function of the normal distríbution.
21

---

Published in Transations on Machine Learning Research (07/2023
Proof of theorem We can prove theorem Since the variance of equation solely depends on the number of time steps, which is equal for all experiments, distributions can only be different because of their means. We start with experiments that are designed according to equation 2a In this case, for distributions to be equal, and, thus, for variables to be non-causal, we need
es (') (A
1-1 Σ (Bultt - = cs (A'"((0))) 1=0
+ -e())))(((
1-1 ABu(tt -
+
==
where 6 E R the unit vector (i.c., a vector with zcros and a single 1 at the ith entry). Since input trajectories are equal, this boils down to
c - -0==((0)=,AA''"00
Essentially this means that the component ij of A needs to be O. This is clearly the case, there is no inftuence of z; on 2j, le., in case variables are non-causal, we have MMD O. The event component ij of A being by chance, even though zy has a causal influence on Fe, has probability 0. Thus, we have that variables non-causal if MMD = 0.
For experiments that are designed according to equation 2b initial states equal and, in case variables are non-causal, we have
Σ 1 Σ A' 1- 1
=
Similar before, we have equal distributions and, thus, MMD 0 if entries in the A'B matrices relating 24 and uj 0, ie., if there is no causal influence. The other direction holds since the event of the relevant entries being 0 by chance has probability 0.
B Further Example for Assumption
We provide a further example to empirically validate assumption 1 In practice, regions of local non-causality are often due to hysteresis effects. For instance, a resting body first needs to overcome static friction. That is, when the velocity is zero, the input signal needs overcome threshold actually influence the velocity. In particular, consider the system
0.01
(kk) ((.1)
( u(k) + u(k) if if + 40V [u(k] 0.1 0.1
(I + 1
0
(17)
(xk)
clse,
where the state 2-dimensional, the input is scalar, and v(k) is a normally distributed random variable with standard deviation 1 10-4 We excile the system with a ramp function starting zero and going until two and learn GP model. For the GP model, we use Gaussian kernel for which we optimize the hyperparameters. We then assume the initial condition of both state variables be zero and compute the MMD of simulated trajectories for 100 step inputs, with u' ranging from 0.05 to 5 and u! -u!. The results confirm our findings from example since also here assumption is satisfied (see figure g
C Further Results of the Robot Experiments
We first provide implementation details for the experiments presented in section D The initial model estimate is obtained by exciting the system with a chirp signal for 30 and using the generated data to learn a linear state-space model (cf. equation 6 with least squares. The obtained matrices are
0.868
0.132
0.754
-0.491
00.075
-0.056
-0.031
0.022
0.132
0.868
0.754
-0.491
0.074
-0.055
-0.031
0.022
Amit
Bunt
-0.132
-0.132
1.754
-0.491
0.074
-0.056
-0.03
0.022
-0.134
-0.134
0.76
0.508
0.075
-0.056
-0.032
0.022
22

---

Published in Transations on Machine Learning Research (07/202)
0.8
0.6
MMD'
0.4
0.2
0.5
1.5
2
2.5
3
3.5
4
1.5
Figure 8: Simulated MMD) for the hysteresis system in equation see that assumption is satisfied.
17
Similar as for erample 7 also here
Initial conditions and input trajectorics for the causality testing experiments are obtained through sophisti- caled guesses, as discussed in section 5.2 The found initial conditions and input. trajectories yield expected need only eight experi- MMDs orders of magnitude above the system's noise level for all joints. Thus,
ments to identify the causal structure. We design input trajectories of 100 time steps for each experiment, repeat the experiment ten times, and use collected data from all experiments for hypothesis testing. While the squared MMD is always positive, the empirical approximation in equation 7 can become negative since it is an unbiased estimate. For the test statistic in equation g we estimate the variance using 100 Monte Carlo simulations and obtain the expected value through a noiseless simulation. We use v 1, but we will see in the results, the empirical MMD is, in all cases, orders of magnitude below or above the threshold. Thus, more consorvative choices of would yield the same outcome.
In tables D and 2 present the results of all causality lesting experiments conducted on the robotic plat form shown in figure 3 As for the results discussed in section O we always have a clear decision on whether to accept or reject the null hypothosis: The MMD found in experiments is always ordors magnitude larger or smaller than the test statistic. Also here, we find that all joints can be moved independently of each other and are affected by exactly one input. When exploiting the revealed causal structure for identifying the system matrices, we obtain
0
0.013
0
0
1
0
0
0
0.007
0
0
Acaus
Beaus
2
0
1
0
0
0.01
0
0
0
1
0
0
0.01
Further Results of the Quadruple Tank Experiments
We discrctize the quadruple tank system with a time-step of 100 ms. We choose A, = 50 em2 for all tanks, 01,2 = 0.242 cm², cm² 03.4 = 0. 0.242 cm², and the valve parameters Ç1.2 0.5. For the initial model learning, we excite the system for 5000 time steps. During excitation, the input is sampled from univariate distribution with an interval [0, 60). We use the generated data to learn GP with Gaussian kernel for each state variable. Having identified model, we follow algorithm [ to identify the causal structure. As for the robot experiments, repeat each experiment 10 times and use 50 simulations to est the standard deviation. We choose 10 in equation to avoid false positives. The results are collected in tables 3 and
23

---

Published in Transations on Machine Learning Research (07/2022)
Table 1: Results of the causal structure identification for robot arm. Causal inftuences of joints on each other.
Joint +
Experimental MMD
Test statistic
Joint
21 + 21
0
1.65 x 10
21 + 22
0
1.79 x 10-4
21 + 23
0
2.39 x 10-4
21 + 2/
2.43 x 10-13
1.61 10
22 + 21
0
5.6 10-7
22 22
2.8 10 18
4.58 10-7
22
0
3.56 x 10-7
22 +
1.82 x 10-13
6.91 X 10-7
33 + 21
0
5.81 X 10
23 + 22
0
4.54 X 10-7
T3 + 23
1.38
10-18
6.16 10-7
33 + JA
10-16 1.29
5.2 10-7
34 21
0
4.99 10-7
24
0
4.66 10
34
9.63 X 10
5.8 X 10-7
34 + 24
-5.44 x 10-15
5.8 x 10-7
Table 2: Results of the causal structure identification for robot arm. Causal influences of inpuls joints.
Input + Joint
Experimental MMD
Test
statistic
01 11
0.04
1.18 10
+ 32
0
5.38 x 10-7
+ I3
0
6.51 x 10
+ 44
0
4.47 x 10-7
3.09 x 10 10
122 + I1
1/2 >> J72
0.04
1.15 10-3
3.51 10-14
5.95 10-7
U2 + .F3
1/2 - 2F.
3.51 10 10
4.67 10
U3 + F1
2.31 X 10
5.26 10
- 422
1.4 x 10
4.28 x 10
03 + 403
0.04
1.18 x 10
d3 - XA
2.25 X 10 12
4.36 x 10
u1 - X1
4.47 x 10 1
5.11 10
U4 - F2
5.69 10-7
u4 - 373
UA + AA
0.04
6.58 x 10
24

---

Published in Transations on Machine Learning Research (07 (07/2022)
Table 3: Results of the causal structure identification for quadruple tank system. Causal influences of the tanks on each other.
Tank + Tank
Experimental MMD
Test statistic
Ø1 + x1
2.78 x 10
9.81 x 10 40-4
21 - 82
7.39 10-6
9.55 x 10-5
31 + 33
1.47 10-5
9.78 10-3
9.45 10-5
1.01 x 10-3
J.4
42 + 21
3.57 10-6
2.35 X 10
+ I2
10-110-1
4.58 X 10
42 + 323
x 10-6
4.31 x 10
42 + J4
3 10-5
1.41 X 10
43 7 I1
1.62 10
1.29 X 10
x3 + X2
5.36 x 10
3.28 10
33 33
10-1 3.44 10-1
2.75 x 10
33 + AAA
9.59 10
3.46 x 10
-
4.28 10
4.04 x 10
+ X2
1.45 x 10-1
8.74 x 10-5
F4 - D3
2.35 x 10-6
6.56 x 10
44 7 JF4
3.57 x 10
1.42 x 10
Table 4: Results of the causal structure identification for quadruple tank system. Causal influences of inputs tanks.
Input +
Experimental MMD
Test statistic
Tank
U1 + 01
1.8 x 10-1
1.01 10
U1 - 52
8.04 10-3
8.74 10-3
+ 43
2.46 x 10-8
3.8 10-5
+ IA
1.56 x 10-1
8.61 x 10
1/2 - 41
7.55 x 10-3
9.86 x 10-5
122 + D22
1.78 x 10-1
1.12 x 10
02 - 373
1.58 x 10-1
3.77 10-5
42 + AAA
5.44 1
4.46 10 10-5
25


Reviews:
[{"review_id": "1iKbpoOUTc", "reviewer": "Reviewer_7Puu", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The proposed algorithm (formulations: 5a, 5b) is an elegant application of the MMD methods and the authors do a good job in explaining how to make this approach scalable. What I found very important is the discussion how to achieve desired states for the presented algorithm by controlling the system which is crucial from the practical point of view (Section 4.2 on controllability). Presented empirical results show that the method is very effective in practice do determine causal relationships.\n\nThe weakness of the approach is that even though elegant, presented algorithm is pretty straightforward given used definition of non-causality combined with Assumption 1 (MMD formulation is almost tautologically implied by Definition 1 and Assumption 1). One of my concerns is that it is not clear at all which dynamical systems f satisfy Assumption 1 or how to test this efficiently. Given the questions regarding how constraining Assumption 1 is, it is not clear to me how widely applicable the methods presented by the Authors are.\n\nThe paper would also benefit from adding few more difficult and high-dimensional dynamical systems setting, in particular to test the scalability of the algorithm and its data efficiency.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The paper presents an algorithm for identifying causal structure of the dynamical systems, predicting which state/control dimensions are causally related and which are not. Determining causality is increasingly important in controlling dynamical systems as providing powerful inductive bias that would be hard to determine from data automatically via deep ML. The algorithm presented by the authors is based on the non-causality definition related to the one from Pearl 1995. The algorithm comes in two different variants: (a) the one where initial conditions vary whereas controls are thee same and (b) the one where initial conditions vary, but the controls differ. The method proposed by the authors uses maximum mean discrepancy (MMD) based metrics (measuring discrepancy between trajectories) to identify causal relationships and discuss how to apply them in the control setting, where states can be achieved only by applying a sequence of controls. In order to scale up the approach, the authors leverage kernel-based unbiased estimator of the squared MMD wit the use of the characteristic kernels. Presented theoretical results are accompanied with empirical evaluation on: (a) the kinematic control of the robot, (b) the quadruple tank process (to test the setting with highly nonlinear dynamics) and some synthetic tasks. ", "claims_and_evidence_yn": "", "claims_explanation": "", "audience_yn": "", "audience_explanation": "", "requested_changes": "Strengthening the experimental section would be very welcome. The authors can also work on improving the presentation of the paper. For instance: (a) the Authors refer to \"do-\" notation in Definition 1, without properly introducing it earlier, (b) while providing an introduction to MMDs the authors refer to characteristic kernel as useful while working with the empirical MMD, but do not provide their definition (just Gaussian kernel as an example). Most importantly, the Authors should comment on the Assumption 1 which is critical for the correctness of the algorithm (the proof that the proposed method is almost immediate if we use it).", "broader_impact_concerns": "I do not see any concerns regarding ethical aspects of this work.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "H60qXEbD9e", "reviewer": "Reviewer_VmF2", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The paper addresses an important problem on causality in the context of control systems. As the paper discusses, most causal identification algorithms are either entirely based on observational data or interventions where one has full control over the variables of interest. Therefore, I believe that an approach that takes controllability restrictions into account is interesting for the machine learning community working on causality. However, the paper also comes with a number of issues which make me believe it's not ready yet for publication.\n\n#### Major issues:\n- 1) **Insufficient theoretical background.** The proposed method is based on Judea Pearl's do-calculus, but no (minimal) background on do-calculus is provided. Although one could argue that is possible to understand most of the paper without this background, the problem setting in Sec. 3.1 gives the impression that it is actually a somewhat necessary background.\n\n- 2) **Notation issues.** The notation in the paper can be quite confusing sometimes, as some of the symbols have their definition coming late or are left implicit/unexplained. For instance, the subscripts and superscripts in Definition 1 have no prior definition. Later on, it can be *assumed* that the subscript $i$ on a state variable $x_i$ denotes the $i$th coordinate of that variable, or perhaps an individual vector-valued variable within a group of vector-valued variables composing the state of the system (the paper does not make it clear). Nowhere in the paper I have seen a formal specification for the control and state spaces. They're only left to be assumed as some possibly arbitrary vector space of some dimensionality. It is also not clear what the superscripts \"I\" and \"II\" represent in Definition 1 in a general context. Furthermore, $\\mathcal{X}_{nc}$ and the corresponding $\\mathcal{U}$ are only loosely explained, but not formally defined.\n\n- 3) **Confusing theoretical statements.** Assumption 1 sounds more like a theorem than an assumption due to the last sentence in its statement. It might need to be split into assumption and theorem/lemma/remark statements, depending on whether it requires a proof or the proof happens to be trivial.\n\n- 4) **Hyper-parameters.** The proposed method has a few hyper-parameters, such as the MMD thresholds $\\delta_1$ and $\\delta_2$, but the paper apparently does not discuss how these hyper-parameters can be tuned or set accordingly in practice for different applications.\n\n- 5) **Lack of experimental baselines.** The experiments section does not provide comparisons against other suitable algorithmic baselines, other than  SINDy in the last synthetic experiment. Although existing work is mostly based on observational data or full interventions, it would be interesting to see how the proposed algorithm fares against other methods in settings where they can be applied. As far as I understand, at least the methods based on purely observational data can be applied to the same problems explored in the paper. An interesting question would then be whether they can identify the same causal relationships.\n\nMinor issues:\n- I believe there's a typo in the definition of $m_2$ in Lemma 2.\n- Kinematic model in experiment 1. Torque is proportional to angular acceleration, not speed, as in the equation $\\dot{\\phi}_i = \\tau_i$, if $\\phi_i$ indeed represents angular position/angle.\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper presents an experimental design method for causal structure identification in general, possibly non-linear dynamical systems which are not necessarily fully controllable. The method is based on comparing the trajectory distributions under different settings for initial conditions and control signals. Theoretical results are derived relating the maximum mean discrepancy (MMD) of the trajectory distributions and causal relations within controllability limitations. A set of experiments demonstrating the method on a physical robot and two synthetic problems is also presented.", "claims_and_evidence_yn": "", "claims_explanation": "", "audience_yn": "", "audience_explanation": "", "requested_changes": "Considering the issues presented above, I suggest:\n- 1. Adding a brief theoretical background on do-calculus in the context of causal identification\n- 2. Revising the notation in the main theoretical results\n- 3. Revising theoretical statements, making sure that notation is properly explained and formally defined\n- 4. Adding a discussion on hyper-parameters tuning\n- 5. Ablation experiments showing how the method behaves as a function of hyper-parameter settings\n- 6. Adding further experimental baselines to experiments in Sec. 6.", "broader_impact_concerns": "No major broader impact concern, other than perhaps adding a brief discussion on the possible consequences of inferring causality in experiments outside of robotic settings, which may involve human agents, for example, case that's a suitable and relevant application.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "l36O5v7hRM", "reviewer": "Reviewer_QpEU", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths. This is a work of possible interest to the community. Writing and presentation are clear and, for the most part, concepts and statements are clearly introduced and treated in a rigorous manner. The idea introduced in the paper is simply yet well founded.\n\nWeaknesses. \n\n(i) I think that knowledge of do-calculus might not be that common among readers. Doing a soft introduction to it might be beneficial to the paper. The do-operator appears in Definition 1 without much introduction (even though it is explained afterwards). Maybe a \"Preliminaries\" section would help the reader more unfamiliar with do-calculus.\n\n(ii) The related work section and the overall framing of this work with respect of the state of the art is somewhat odd. Most of the works cited, while important works, are quite old. While reading the manuscript I had a hard time placing this work with respect to more recent works addressing the problem of identifying causality. Overall, it feels it feels like a work disconnected from recent work. I'd suggest a better reframing of the Related Work section to put better emphasis on how this work relates to recent state of the art.\n\n(iii) Why is the maximum mean discrepancy used as the measure of similarity (instead of something like the KL divergence). This is not clearly explained/motivated. A somewhat vague justification is given for the case of reproducing kernel Hilbert spaces, but I'd suggest clarifying this further.\n\n(iv) Overall, the information presented in the numerical results is simplistic. For example, a motivation for learning causality is to reduce model complexity, but this is not really evaluated in the numerical results. In the same sense, for the quadruple tank process, only the causality tests are shown (and extended in the Appendix). However, it would be interesting to see how a reduced model learned from the causality learning process behaves. I'd suggest extending the numerical results with more insightful figures and discussion.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "In this paper, the authors study the problem of identifying causality in dynamical systems. To this end, they introduce an algorithm based on measuring maximum mean discrepancy. Implicitly, the method considers controllability, making it a viable method to identify causality in control system. The method is evaluated in simulations consisting of a robotic arm and a quadruple tank process.", "claims_and_evidence_yn": "", "claims_explanation": "", "audience_yn": "", "audience_explanation": "", "requested_changes": "Addressing the previously mentioned weaknesses should be sufficient.", "broader_impact_concerns": "None.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: Learning the causal structure of an unknown dynamical system and exploiting it for data-efficient control is a very well motivated and important problem, with many applications of cross-disciplinary interest. This paper develops a class of MMD based methods to compare trajectory distributions under varying initial conditions/control sequences to infer causality, which implicitly accounts for controllability of the underlying system -- the reviewers agree this is an elegant approach overall, presented clearly and rigorously in the paper. Several concerns and issues were raised around comprehensiveness of the experimental results, but they were diligently addressed during the review period and the AE believes the expansions described in the author responses are sufficient. For final revision, all reviewers requested a more complete and self-contained introduction to "do-calculus", related work that better covers more recent works, and a stronger commentary around what classes of nonlinear systems may or may not be captured by Assumption 1. 
Decision: Accept (poster)