<figure><img src="image_1.png" /></figure>

HAL
open science

# Dynamics-Based Identification of Hybrid Systems using Symbolic Regression

Swantje Plambeck, Maximilian Schmidt, Goerschwin Fey, Audine Subias, Louise Travé-Massuyès

► To cite this version:

Swantje Plambeck, Maximilian Schmidt, Goerschwin Fey, Audine Subias, Louise Travé-Massuyès. Dynamics-Based Identification of Hybrid Systems using Symbolic Regression. DSD/SEAA Conference, Aug 2024, Paris, France. hal-04794430

HAL Id: hal-04794430

https://hal.science/hal-04794430v1

Submitted on 20 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
# Dynamics-Based Identification of Hybrid Systems using Symbolic Regression

Swantje Plambeck, Maximilian Schmidt, Goerschwin Fey
Institute of Embedded Systems

Hamburg University of Technology
Hamburg, Germany

firstname.lastname@tuhh.de

Audine Subias*, Louise Travé-Massuyès†

*† LAAS CNRS, * INSA

Université de Toulouse
Toulouse, France

*subias@laas.fr, †louise@laas.fr

*Abstract*—Symbolic regression has shown potential in the identification of physical systems. Hybrid systems, which combine both continuous and discrete behavior, are a relevant extension of purely physical systems, used in many fields, including robotics, biological systems, and control systems. However, due to their complexity, finding an accurate model is a challenge. This paper presents a novel approach to learning models of hybrid systems using symbolic regression. Our method leverages the power of genetic programming to automatically discover accurate and interpretable mathematical models in the form of hybrid systems from observed data. Symbolic regression detects transitions between different continuous behavior of a system directly based on the dynamics, instead of pure distances of observed trajectories. Furthermore, models generated by symbolic regression can be used to predict future system behavior, detect anomalies, and identify the underlying dynamics of the system while providing a human-readable representation. Our results demonstrate that symbolic regression can effectively identify the underlying dynamics of a real system represented in a hybrid model, providing a valuable tool for system identification and diagnosis.

*Index Terms*—Hybrid Systems, Symbolic Regression, System Identification

## I. INTRODUCTION

Hybrid systems are abstract models of systems that exhibit both continuous and discrete behavior. They are used to model a wide range of systems, including cyber-physical systems, manufacturing systems, and many more [1, 2]. Due to their inherent combination of continuous and discrete behavior, the identification of hybrid systems is a challenging task. Nevertheless, accurate abstract models are essential for control, simulation, and verification of these systems. Hybrid systems have a finite number of modes representing a specific inherent dynamic of the system. Within a mode, the continuous system dynamics are defined by flow functions over the system variables. An event occurs based on guard conditions on the system variables. The event triggers a transition from one mode to the next.

We propose a novel approach to detect the transitions between modes and the identification of the continuous dynamics themselves using Symbolic Regression (SR). Existing hybrid system identification methods including previous work of the authors introduce a general procedure as shown in Figure 1 that identifies modes and their flow functions to create a model [3]. Our approach makes contributions in

<figure><img src="image_2.png" /><figcaption>Fig. 1: Steps of Hybrid System Identification</figcaption></figure>

improving the steps (1) to (3) for detection of transitions, grouping of observations with same dynamics into modes, and identification of modes via symbolic expressions. SR is a powerful tool for the identification of continuous dynamics, where it has been shown to be able to identify complex behavior from data [4, 5, 6]. Within the identification of hybrid systems, we use this ability to detect not only the expressions describing dynamic behavior, but also to 1) separate data segments that cannot be represented by the same expression, as well as 2) group data segments that are represented by the same expressions. Finally, SR provides human-readable mathematical expressions representing the system dynamics, which is essential for the interpretation of the identified models and supports diagnosis and fault explanation for the system.

Our approach aims to detect mode transitions directly based on the dynamics of the observed data. This goes beyond existing identification strategies, which focus on similarities of observations to detect whether observations belong to the same or different modes of a hybrid system [7, 8, 9]. Our approach assumes that one discrete mode is governed by a single symbolic expression and is able to handle cases where multiple modes are necessary to represent the system dynamics.

The main contributions of this paper are:

* A novel approach for the detection of transitions between different modes of a hybrid system using SR, separating observed trajectories into segments,

* a grouping of segments to identify the different modes of a hybrid system re-using SR,

* an evaluation of the proposed approach on example systems including a physical simulation of a two-state power converter and a two tank system.

The first two contributions are the main focus of this paper, giving an algorithmic approach for the identification of hybrid systems using SR. We start from a small window of the
observed data, learn a symbolic expression, and then iteratively
enlarge this window until a transition is detected. A transition
is detected as soon as the fitness of the learned expression
drops. The segments between detected transition points are
then grouped to identify the modes of the hybrid system. For
this grouping step, we re-use SR to learn expressions on unions
of the previously detected segments. Whenever the loss of
the combined segments declines compared to the individual
segments, the segments are grouped. The expression that is
found on this segment identifies the mode.

The paper is structured as follows: in Section II, we review
related work on system identification, specifically for hybrid
systems and symbolic regression. In Section III, we introduce
the necessary formal definitions. In Section IV, we present our
approach for the identification of hybrid systems using SR –
covering the detection of transitions, grouping of segments,
and mode identification. In Section V, we evaluate our ap-
proach on two representative examples. Finally, in Section VI,
we conclude the paper.

II. RELATED WORK

A. *Symbolic Regression for System Identification*

Symbolic Regression (SR) is a method for regression and
system identification that aims to find a symbolic expression
that matches a given data set. Contrary to traditional regression
methods, SR is not restricted to a specific set of functions or
structure of an expression like, e.g., polynomial regression, but
uses a set of basic operators to construct complex expressions
freely. The development of SR has been supported by the
advancement of genetic programming, which is commonly
used to implement SR algorithms [5]. In addition to genetic
programming, there exist further approaches for SR, e.g., using
deep reinforcement learning [10] or lattices [11].

SR has been applied to a wide range of problems, including
the identification of physical concepts from data [4], hybrid
dynamical systems [12], for general system identification
problems [6], mining the expression of diagnosis indicators
[13], and even for the identification of a machine learning
algorithm [14]. Thus, SR is one of several methods for
system identification together with a broad range of other
identification methods such as linear or nonlinear regression,
neural networks, or kernel-based methods [15]. Schoukens
and Ljung [16] provide a comprehensive overview on these
methods.

B. *Identification of Hybrid Systems*

In this paper, we focus on the identification of hybrid
systems, which are often equivalently considered as hybrid
automata. There exist several approaches on the identification
of hybrid systems, e.g., using clustering methods [7], involving
machine learning methods with neural networks [9], or linear
matrix inequalities [8]. All of these methods follow a similar
procedure consisting of steps for the detection of transitions,
grouping of segments, and mode identification. For the detec-
tion of transitions, most methods use a similarity measure on
the observations, e.g., based on distance [7], in the frequency

domain [9], or based on windows of the observations [8].
These approaches are not aware that the same dynamics
can generate complex trajectories whose observations lie in
different value ranges. In contrast, our approach focuses on
the detection of transitions directly based on the dynamics
appearing in the observed data.

### III. PRELIMINARIES

In this section, we introduce the basics of Symbolic Re-
gression (SR) and formally define hybrid systems and system
observations.

A. *Symbolic Regression*

SR is a machine learning technique that aims to find a
mathematical expression describing the relationship between
input and output variables. The goal is to find a function r
that maps input variables i to an output variable o such that
o = r(i). The function r is represented as a mathematical
expression in terms of basic functions and operators [17].
The search for the best symbolic expression is guided by a
fitness metric that evaluates the quality of candidate symbolic
expressions on a learning data set.

Learning algorithms typically represent the expression as a
tree structure, where the nodes represent the operators and the
leaves represent variables or constants. The search space 𝒬 is
the set of possible tree structures defined by a set of operators,
as well as the set of variables and constants. The search space
can be constrained by limiting the maximum depth of the tree,
the number of nodes, or the set of operators [10, 18].

In the scope of this work, we use the framework PySR for
SR, which uses Genetic Programming (GP) [17, 18]. Using
GP, the search for the best symbolic expression is performed
by evolving a population of candidate expressions over so-
called *generations*. The population is initialized with random
expressions. In each generation, the population is evaluated
using the fitness function. New solutions are generated by
applying genetic operators, such as mutation and crossover, to
the best candidates. The process is repeated for a predefined
number of generations.

The convergence speed of GP has a high variation and
convergence may take long, especially for complex problems.
There are no formal guarantees that the algorithm will find
the optimal solution. Conditions for the convergence of GP
are still an active field of research [5].

Further challenges for model learning – that are not specific
to GP – rise from insufficient data quality. Especially on
noisy data, the actual underlying function might not be found.
However, GP is known to find good approximations of the
actual function [4, 5]. Still, the amount of data and the
range of captured behavior are crucial for the success of
the algorithm. Additionally, the complexity of the searched
expression influences the difficulty of the learning task. Here,
we assume that the data is generated by a real system that can
be abstracted by a hybrid system. Furthermore, the collected
data is sufficient to capture the flow functions of the modes
in the hybrid system abstraction.
## B. Hybrid Systems

**Definition 1 (Hybrid System [19]):** A Hybrid System is defined by a 5-tuple $(X, Q, \mathcal{F}, \mathcal{T}, \Sigma)$ where

* $X = \{x_1, x_2, ..., x_n\} = I \cup O \cup S$ is the set of system variables which consist of input variables $I$, output variables $O$ and state variables $S$. Derivatives of variables may form individual variables in $X$.

* $Q$ is the set of modes.

* $\mathcal{F}$ is the set of flow functions. Every flow function $f_q \in \mathcal{F}$ defines the change of the state variables $S$ as well as the current output variables $O$ over the continuous time $t$ within the mode $q$ based on the current values of state variables and the external inputs $I$, i.e., $[O(t), \dot{S}(t)] = f_q(S(t), I(t), t)$.

* $\mathcal{T}: Q \times \Sigma \to Q$ defines transitions between modes $Q$. A transition is triggered if the corresponding event $\sigma \in \Sigma$ is active.

* $\Sigma$ is a set of events leading to transitions between modes. Each event is guarded by a set of conditions on the variables in $X$. A transition is triggered if the conditions are met.

A hybrid system, according to this definition, combines discrete and continuous behavior. Discrete behavior is captured by the discrete modes $Q$ and transitions $\mathcal{T}$, while the flow functions $\mathcal{F}$ describe the continuous dynamics. In contrast to finite automata, we exclude instantaneous state changes with hybrid systems, but consider a positive residence time in each of the hybrid modes. Transitions are usually triggered by conditions on the variables or discrete control signals. In the scope of this paper, we focus on the identification of modes and flow functions while excluding the construction of conditions. Within this scope, dynamics define the physical behavior of the real system. In the abstraction given by the hybrid model, the dynamics are represented by the flow functions of the modes. Observations of dynamics, i.e., changes in the values of the variables, are considered as trajectories which are introduced in the following.

## C. Observations

We define a trajectory as an observation in form of a multidimensional time series of the variables $X = \{x_0, x_1, ..., x_n\}$. The trajectory is a sequence of samples of $X$ at distinct sampling points in continuous time.

The input to the identification algorithm is a set of trajectories. Each trajectory covers multiple changes between the modes and stems from observing the real system for a finite time interval.

# IV. IDENTIFICATION OF HYBRID MODES

## A. Overview & System Properties

**Definition 1** specifies hybrid systems in a general form. The most characteristic property of hybrid systems is the combination of continuous dynamics (defined by the flow functions $\mathcal{F}$) with discrete modes ($Q$). To learn both parts, identification of hybrid systems includes multiple subproblems

as shown in Figure 1 and as similarly introduced in Saberi et al. [20]:

1. detection of discrete transitions between dynamics, separating the trajectories into segments,

2. grouping of segments with same dynamics, forming discrete modes,

3. identification of the continuous dynamics for each mode, i.e., the flow functions of $\mathcal{F}$,

4. model construction, i.e., accumulation of the results in a single hybrid system.

In this section, we propose a novel approach for steps 1. to 3. in which steps 2. and 3. are combined. We hence present two algorithms: one algorithm to detect the transitions between modes of a hybrid system using SR and, further on, a second algorithm to group and identify the flow functions of modes with SR using the segmented trajectories. We assume the following properties of the real system and its abstraction as a hybrid system:

1) Each mode can be described by a symbolic expression using a suitable set of basic operators. Within each mode, the variables evolve continuously. Discontinuities are only possible at mode transitions. We, thus, require flow functions to be continuous, i.e.,

$$ \lim_{t \to t_a} f_q(t) = f_q(t_a), \forall t_a \in \mathbb{R}, \forall f_q \in \mathcal{F} \quad (1) $$

where $t$ is the continuous time.

2) The system may have multiple dynamics, i.e., more than one mode is needed to represent the system dynamics ($|Q| > 1$). This introduces the need to detect transitions between modes.

Our approach always considers a single output for all steps of the identification procedure. The extension to multiple outputs is straightforward, but requires additional consideration for the detection of segments, as we have to combine detected transitions for all outputs. Further, we assume completeness and distinctness of training data with respect to the desired model. This is a common assumption, usually ensured by data quality requirements during data aggregation. We assume that no binary or categorical variables which directly indicate mode transitions are available. Furthermore, no discrete variables are involved in the construction of the continuous flow functions.

## B. Detection of Transitions

Our approach for the detection of transitions between different modes of a hybrid system uses SR. Algorithm 1 describes the procedure. Given an observed trajectory of the system, we process over this trajectory using a window (Line 5). Initially, this window covers a fixed initial length at the beginning of the trajectory. As long as the segmentation criterion is fulfilled, i.e., the loss of a learned expression does not significantly drop, the window is extended to the right (Line 8). As soon as the segmentation criterion is not fulfilled anymore, a transition between modes is detected and the window is stored as a segment in the set $\mathcal{T}$. In the inner loop, the symbolic expression is learned incrementally, i.e., $n_{update}$-many iterations of
the SR are performed. This ensures that the knowledge about sub-expressions or constants from past SR generations are re-used on window extension. Still, an adaption of the previous candidate expression is performed in case that the new data from the extended window provides additional information on the dynamical behavior. For example, an initially almost linear change may actually turn out to be part of a higher order expression. Whenever a transition is found, learning restarts (Line 15) and creates a new window starting at the end of the previously identified window. A new expression is learned from scratch.

**Data:** trajectory

**Result:** $\mathcal{T}$, expressions

1. loss ← deque[l<sub>hist</sub>];

2. i<sub>start</sub> ← 0; i<sub>end</sub> ← l<sub>init</sub>; n ← n<sub>init</sub>;

3. **while** i<sub>end</sub> < len(trajectory) **do**

4. &nbsp;&nbsp;&nbsp;&nbsp;**while** loss is empty or segmentCriterion(loss) **do**
5. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;window ←
6. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;trajectory[i<sub>start</sub>, min(i<sub>end</sub>, len(trajectory))];
7. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;exp, fit ← learnExpression(window,n);
8. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;loss.append(fit);
9. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;i<sub>end</sub> ← i<sub>end</sub> + l<sub>step</sub>;
10. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;n ← n<sub>update</sub>;

10. **end**

11. $\mathcal{T} \leftarrow \mathcal{T} \cup \{\text{window}[0, \text{end} - l_{\text{step}}]\}$;

12. expressions.append(exp);

13. loss.clear();

14. i<sub>start</sub> ← i<sub>end</sub> − l<sub>step</sub>; i<sub>end</sub> ← i<sub>start</sub> + l<sub>init</sub>;
15. n ← n<sub>init</sub>;

15. resetSR;

16. **end**

**Algorithm 1:** Detection of Mode Transitions

The algorithm is parametrized by:

* $l_{init}$: the initial window size when learning an expression from scratch,

* $l_{step}$: the step-width for extending the window in the inner loop,

* $n_{init}$: the number of iterations of SR when learning an expression from scratch,

* $n_{update}$: the number of iterations of SR for updating the learned expression on an extended window,

* $l_{hist}$: the history of loss considered for the segmentation criterion.

Furthermore, SR uses additional parameters for learning an expression, as introduced in Section III-A. Finally, the functions `segmentCriterion` and the metric used to rate the deviation between the learned expression and the original data are parameters of the algorithm. Here, we use the mean-square error loss to evaluate the deviation. We apply a thresholded depression-condition as follows:

$$ (\text{loss}[\text{end}] < \tau) \text{ or } (\text{loss}[\text{end}] \le \text{loss}[\text{end} - 1]), \quad (2) $$

where $\tau$ is a stagnation-threshold on the latest loss-value. The stagnation-threshold is chosen based on the required or

assumed loss of learned flow functions.

### C. Grouping and Identification of Modes

Segments whose dynamics follow the same symbolic expression refer to the same mode of the system. Thus, the next step involves grouping segments with identical dynamics. Algorithm 2 shows the pseudocode for segment grouping using SR. The input to the grouping is the set of detected segments $\mathcal{T}$ from the previous step. The first segment is a first candidate group, as stated in Line 1. Afterward, we iterate for every segment in $\mathcal{T}$ in Line 2 over all known groups in Line 4. The current segment and the current group are combined to one data set and an expression is learned with SR (see Line 5). Re-learning the expressions during grouping is computationally expensive, but allows adapting to minor variations in the dynamics. The number of generations $n_{group}$ for SR is a parameter of the algorithm. If the loss of the learned expression is small, the current segment is included in the current group (see Lines 6 and 7). If no matching group is found for the current segment, the segment forms a new group, as stated in Line 13.

**Data:** $\mathcal{T}$

**Result:** $G$, expressions

1. $G \leftarrow \{S[0]\}$;

2. **for** $s \in \mathcal{T}$ **do**

3. &nbsp;&nbsp;&nbsp;&nbsp;groupFound ← False;

4. &nbsp;&nbsp;&nbsp;&nbsp;**for** $g \in G$ **do**

5. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;exp, fit ← learnExpression($s \cup g, n_{group}$);

6. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**if** groupingCriterion(fit, s, g) **then**

7. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;g ← $s \cup g$;
8. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;expressions[$g$] ← exp;
9. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;groupFound ← True;
10. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;break;

11. &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**end**

12. **end**

13. **if** not groupFound **then**

14. &nbsp;&nbsp;&nbsp;&nbsp;G.append(s);

15. **end**

16. **end**

**Algorithm 2:** Grouping of Modes

As for the segmentation, the loss is used as the metric for the deviation between learned expression and original data. The `groupingCriterion` is a parametrization opportunity for the algorithm. We use the following grouping criterion:

$$ \text{loss}(s \cup g) \le \varphi \cdot \text{loss}(g), \quad (3) $$

where $\text{loss}(g)$ gives the loss of the group observed before combining with the segment $s$ and $\text{loss}(s \cup g)$ is the loss after combination – the smaller, the better. The parameter $\varphi$ is a relaxation parameter.
## V. EXPERIMENTS

In this section, we present empirical results on example systems. In the discussion of the performance of our approach, we consider the detection accuracy for transitions, compare learned expressions and ground truth flow functions, and evaluate the overall accuracy in predicting the trajectory of a system.

### A. Example Systems

We present empirical results on two example systems. The first example is a well-known benchmark system, and the second example is a simulation of a real-world system. The systems are described in the following.

1) **Two Tank System:** The two tank system [21] is a benchmark system for hybrid systems. It consists of two tanks, with a pump pumping water into the first tank and a valve $V_b$ that can be opened to let water flow from the first to the second tank. The system has a controller to regulate the height of the water in the first tank. The following differential equation describes the system:

$$ \dot{h}_1 = \begin{cases} \frac{Q_p - C_{vb} \cdot \sqrt{h_1 - h_2}}{A}, & \text{if } h_1 > h_2, V_b \text{ open} \\ \frac{Q_p + C_{vb} \cdot \sqrt{h_2 - h_1}}{A}, & \text{if } h_1 \le h_2, V_b \text{ open} \\ \frac{Q_p}{A}, & \text{if } V_b \text{ closed,} \end{cases} \quad (4) $$

where $h_1$ and $h_2$ are the heights of the water in the first and second tank, respectively, $\dot{h}_1$ is the derivative, i.e., the change in the water level, of the first tank. $Q_p$ is the flow rate of the pump, $C_{vb}$ is the valve conductance, and $A$ is the cross-sectional area of the first tank. The simulation involves noise and a realistic controller which applies a zero-order hold to the observed height of the tank. The noise is additive noise with a maximum amplitude of $10^{-6}$ for all sensors.

For SR to learn the two tank example, we use addition, subtraction, multiplication, and division as well as the square-root as basic operators. The inflow $Q_p$ and the height of the two tanks $h_1$ and $h_2$ as well as the constants $C_b$ and $A$, and the pre-calculated variable $h_{sqr} = \sqrt{|h_1 - h_2|}$ are given as variables to the learner to learn the derivative $\dot{h}_1$.

2) **Power Converter:** The power converter [22] is a real-world system that has a controlled input voltage $v_s = \sigma V_g$ where $V_g$ is a constant input and $\sigma$ is a switching variable which can be either 1 or -1. In addition to the voltage source, the circuit consists of a capacitor $C$, an inductance $L$, and a resistor $R$. The circuit can be either a parallel or a series configuration, as shown in Figure 2. As presented in [23], we use a transformed coordinate system to describe both configurations with the same equations. The following differential equation describes the system:

$$ \dot{w} = \begin{pmatrix} 0 & \alpha \\ -\alpha & -\beta \end{pmatrix} w + \begin{pmatrix} 0 \\ \alpha \end{pmatrix} \sigma, \quad (5) $$

where $\alpha = (\sqrt{LC})^{-1}$ and $\beta = (RC)^{-1}$ (parallel case) or $\beta = R \cdot L^{-1}$ (serial case). The transformed coordinates $w = [w_1, w_2]^T$ are constructed from the quantities in Figure 2 as

$$ w_1 = \frac{v_c}{V_g} \quad \text{and} \quad w_2 = \frac{1}{V_g} \sqrt{\frac{L}{C}} \cdot i_c. $$

Our goal is to model the state variable $\dot{w}_2$. The inductance, capacity, and resistor have the values $R = 400 \, \Omega$, $L = 8 \, \mu H$, $C = 10.5 \, \text{nF}$, and $V_g = 20 \, \text{V}$.

For SR, the power converter uses addition, subtraction, multiplication, and division as well as the square-root as basic operators. Variables for learning the derivative of $w_2$ are $w_1$, $w_2$, and the continuous time $t$ of a sampling point. Constants are not given as variables, but are estimated by the learner. The to be learned derivative $\dot{w}_2$ is numerically calculated from the data. Different from the two tank example, the simulated data does not involve noise.

3) **General Setup:** As already stated, for the first example, we provide constants as variables to the symbolic regression learner such that the value of constants is not estimated during learning. By this, we focus on the analysis of the algorithmic approach for identification of hybrid systems instead of analyzing the capabilities of SR to estimate constants. Furthermore, for the two tank system, we also provide the pre-calculated term $\sqrt{|h_1 - h_2|}$ as a variable to the learner. This is done to reduce the complexity of the expression that has to be learned. Practically, providing such pre-calculated terms can be interpreted as partial knowledge about the physical system and the operators and terms that are relevant for its behavior.

For both examples, we use the thresholded segmentation criterion from Equation (2) and the mean-square error loss as the loss metric. The experiments show that, using Equation (2), the threshold is often more significant in detecting the transitions than the decline in accuracy. Thus, transitions are mainly detected when the loss is above the specified threshold. The grouping criterion for both examples is the one from Equation (3). The parameters for the learning algorithm for the different examples are shown in Table I. The results provide runtime information for the segmentation and grouping steps. For these measurements, experiments are done on the same personal computer (Intel i7-9750H 2.6GHz, 16 GB RAM, charging).

### B. Detection of Transitions

In a first step, we analyze the transition detection properties of our approach. For this, we compare the detected transitions with the ground truth transitions. Table II shows the results for the two example systems. Transitions are considered correct if they lie in a distance of max. 10 sampling points to the ground truth transition. The columns show the number of detected transitions, the number of true positive transitions, and the number of false positive transitions. The true positive detection shows the accuracy of the detected transitions. The false positive detection indicates the precision. The results show that the approach detects the transitions in both systems with high precision. Nearly all transitions are detected perfectly.
<table><caption>TABLE I: Example Systems</caption><thead><tr><th>System</th><th><i>l</i><sub>init</sub></th><th><i>l</i><sub>step</sub></th><th><i>n</i><sub>init</sub></th><th><i>n</i><sub>update</sub></th><th><i>l</i><sub>hist</sub></th><th><i>τ</i></th><th><i>n</i><sub>group</sub></th><th><i>φ</i></th></tr></thead><tbody><tr><td>Power Converter</td><td>100</td><td>20</td><td>100</td><td>20</td><td>1</td><td>1 · 10<sup>-6</sup></td><td>20</td><td>1</td></tr><tr><td>Two Tank System</td><td>50</td><td>10</td><td>40</td><td>5</td><td>1</td><td>1 · 10<sup>-7</sup></td><td>8</td><td>1.6</td></tr></tbody></table>

<figure><img src="image_3.png" /><figcaption>Fig. 2: Series and Parallel Representation of the Power Converter [22]</figcaption></figure>

<table><caption>TABLE II: Detected Transitions</caption><thead><tr><th>System</th><th>Detected Transitions</th><th>TP in %</th><th>FP in %</th><th>Runtime Segmentation</th></tr></thead><tbody><tr><td>Power Converter</td><td>4</td><td>100</td><td>0</td><td>664.9 s</td></tr><tr><td>Two Tank System</td><td>27</td><td>92.9</td><td>3.7</td><td>520.7 s</td></tr></tbody></table>

For the converter, the detected transitions are shown on the left of Figure 3. For the two tank example, some false positive transitions are detected. This is due to the fact that SR is not always able to adapt the expression on the data adequately, when the modeling window is extended. Instead, additional transitions may be detected. For example, a sinus function may be approximated by a linear function in the beginning, but as soon as more data is added, the candidate expression has to be adapted. Currently, the approach is biased by the initial candidate expression and may detect a transition where there is none. This is a limitation of the approach that is addressed in future work by allowing a re-learning from scratch after the detection of a transition.

### C. Identification of Modes

In the next step, we compare the identified groups and expressions with the ground truth expressions. Table III shows the results for the two example systems. The columns show the number of detected modes and the expressions detected for these modes.

For both examples, the expressions are close to the ground truth expressions. Still, some differences are present, e.g., the power converter uses a square root in the expression, while the ground truth expression does not. As the value range for this term is close to one, the influence of the square root is negligible.

Our approach identifies expressions that are close to the ground truth expression. These expressions usually involve the most dominant terms of the ground truth expression, but may not exactly represent the less dominant parts. This is visible for the two tank system, where the approach is able to identify the dominant influence of the pump defined by the term $Q_p/A$, while the influence of water transfer between the two tanks is identified only for the second, but not for the third mode. The first two modes identified for the two tank example represent the behavior accurately, while the third mode is actually not

needed and an artifact of the grouping step. For the other example, the grouping is correct, which can be seen on the right of Figure 3 by the two colors indicating the two groups.

The runtime of both, segmentation and grouping steps, as given in Table II and Table III is moderate. The algorithm scales especially with the number of evolutions of the symbolic regression, while a single evolution scales especially with the amount of learning data. The grouping step is much faster than the segmentation step, because the step width of window extensions is small compared to the average segment length. Thus, many iterations of the internal loop in Algorithm 1 are needed. The parameter $l_{step}$, thus, is a trade-off between runtime and detection accuracy for transitions.

### D. Accuracy of Predicted Trajectories

In the last step, we compare the predicted trajectories with the ground truth trajectories. For this, we present the mean-square error loss of the predicted trajectories. The results are shown in the 4th column of Table III. The column shows the loss of the predicted trajectories compared to the observed ones, where the loss is the mean square error between the predicted and the ground truth trajectory. Even though the previous section shows that the expressions are not exactly identical to the ground truth expressions for some examples, the predicted trajectories are quite accurate.

Figure 4 shows the predicted and ground truth trajectories for the power converter and the two tank system. The predicted trajectories are close to the ground truth trajectories. Especially, the power converter shows a good fit of the predicted trajectories to the ground truth trajectories. The two tank system shows a good fit in most parts, with two interesting exceptions. The first exception is in the first part of the trajectory. While the ground truth trajectory is very smooth, the predicted trajectory has a small jump. This results from a false positive transition that is detected in the beginning of the trajectory, where SR was not able to adapt a learned expression
<figure><img src="image_4.png" /><figcaption>Fig. 3: Results for the Detection of Segments and Grouping for the Power Converter</figcaption></figure>

<table><caption>TABLE III: Identified Modes</caption><thead><tr><th>System</th><th>Groups</th><th>Learned Expressions</th><th>Loss</th><th>Runtime Grouping</th></tr></thead><tbody><tr><td rowspan="2">Power Converter</td><td rowspan="2">2</td><td>Group 1 – 3.37 · 10<sup>-3</sup> – 4.71 · 10<sup>-3</sup> · w1</td><td rowspan="2">4.76 · 10<sup>-7</sup></td><td rowspan="2">45.8 s</td></tr><tr><td>Group 2 – -6.28 · 10<sup>-3</sup> · √w1 + 0.396</td></tr><tr><td rowspan="3">Two Tank System</td><td rowspan="3">3</td><td>Group 1 – <sup>Q<sub>p</sub>/A</sup>&frasl;<sub>h<sub>1</sub> - h<sub>2</sub></sub></td><td rowspan="3">5.26 · 10<sup>-6</sup></td><td rowspan="3">163.6 s</td></tr><tr><td>Group 2 – (-C<sub>vb</sub> · √|h<sub>1</sub> - h<sub>2</sub>| + Q<sub>p</sub>)/A</td></tr><tr><td>Group 3 – (-C<sub>vb</sub> + Q<sub>p</sub>/h<sub>1</sub>)/A</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>Fig. 4: Predicted and Ground Truth Trajectories</figcaption></figure>

to newly seen data in the transition detection step. The second exception is that some transitions seem to be missed. In fact, all the transitions are detected, but the segments are not grouped correctly. Especially for small value ranges, the fluctuation on calculated losses may lead to a wrong grouping of segments. This could be overcome by grouping always the group with minimal loss, instead of grouping with the first group that has a loss below a given threshold. This comes with a higher computational effort as SR has to be performed for all possible groupings.

For both examples, the approach identifies the structure of the hybrid systems from data with high accuracy. For the power converter the trajectories match the original observations especially well, while for the two tank example some deviations result from the incorrect grouping. Future improvements of the algorithmic procedure may overcome

the limitations and allow for a more flexible and robust identification of hybrid systems from data.

## VI. CONCLUSION

In this paper, we have presented a novel approach to hybrid system identification. We have shown that symbolic regression can be used to identify the structure of hybrid systems from data. We have proposed a new method that uses symbolic regression to iteratively adapt to observed dynamics. Thus, being able to identify similar dynamics on dissimilar observations, as well as transitions between different dynamics on smooth transitions. We have evaluated our method on two example systems: a power converter, and a two tank system. Our results show that our method can accurately identify the structure of hybrid systems from data. Our approach learns a hybrid system without the need for prior knowledge of the
system structure, except for a guess of the involved operators
used for symbolic regression. The output of the identification
process is an explainable and interpretable model which makes
the system understandable to other stakeholders. The model,
thus, serves for prediction of future system behavior and can
be used for control, fault detection, and other applications.

Future work involves advancing the algorithms for better
adjustments on window extension, parameter studies and an
extension of the approach to include the identification of events
and overall model construction.

ACKNOWLEDGMENT

This work was supported by a fellowship of the Ger-
man Academic Exchange Service (DAAD) and the ECIU
Universities. The research is partially funded by the BMBF
project AGenC no. 01IS22047A. It is also supported by
ANITI through the French "Investing for the Future – P3IA"
program under the Grant agreement n°ANR-19-P3IA-0004.
Furthermore, we would like to thank Nicola Zaupa and Luca
Zaccarian from LAAS CNRS for their support and comments
on the power converter example.

REFERENCES

[1] R. Alur, “Principles of Cyber-Physical Systems,” MIT Press, Tech. Rep., 2015.

[2] J. Lunze and F. Lamnabhi-Lagarrigue, *Handbook of Hybrid Systems Control: Theory, Tools, Applications*. Cambridge: Cambridge University Press, 2009.

[3] S. Plambeck, A. Bracht, N. Hranisavljevic, and G. Fey, “Famos– fast model learning for hybrid cyber-physical systems using decision trees,” in *International Conference on Hybrid Systems: Computation and Control*. Association for Computing Machinery, 2024.

[4] M. Schmidt and H. Lipson, “Distilling free-form natural laws from experimental data,” *Science*, vol. 324, no. 5923, pp. 81–85, 2009.

[5] G. Kronberger, L. Kammerer, and M. Kommenda, “Identification of dynamical systems using symbolic regression,” *CoRR*, vol. abs/2107.06131, 2021. [Online]. Available: https://arxiv.org/abs/2107.06131

[6] S. Gaucel, M. Keijzer, E. Lutton, and A. Tonda, “Learning dynamical systems using standard symbolic regression,” in *Genetic Programming*, M. Nicolau, K. Krawiec, M. I. Heywood, M. Castelli, P. García-Sánchez, J. J. Merelo, V. M. Rivas Santos, and K. Sim, Eds. Springer Berlin Heidelberg, 2014, pp. 25–36.

[7] N. Barbosa Roa, L. Travé-Massuyès, and V. H. Grisales-Palacio, “Dyclee: Dynamic clustering for tracking evolving environments,” *Pattern Recognition*, vol. 94, pp. 162–186, 2019.

[8] X. Yang, O. A. Beg, M. Kenigsberg, and T. T. Johnson, “A framework for identification and validation of affine hybrid automata from input-output traces,” *ACM Transactions on Cyber-Physical Systems*, vol. 6, no. 2, pp. 1–24, 2022.

[9] O. Niggemann, B. Stein, A. Vodencarevic, A. Maier, and H. Kleine Büning, “Learning behavior models for hybrid timed systems,” *AAAI Conference on Artificial Intelligence*, vol. 26, no. 1, pp. 1083–1090, 2012.

[10] B. K. Petersen, M. L. Larma, T. N. Mundhenk, C. P. Santiago, S. K. Kim, and J. T. Kim, “Deep symbolic regression: Recovering mathematical expressions from data via risk-seeking policy gradients,” in *International Conference on Learning Representations*, 2021.

[11] K. R. Broløs, M. V. Machado, C. Cave, J. Kasak, V. Stentoft-Hansen, V. G. Batanero, T. Jelen, and C. Wilstrup, “An approach to symbolic regression using feyn,” *CoRR*, vol. abs/2104.05417, 2021. [Online]. Available: https://arxiv.org/abs/2104.05417

[12] D. L. Ly and H. Lipson, “Learning symbolic representations of hybrid dynamical systems,” *Journal of Machine Learning Research*, vol. 13, no. 115, pp. 3585–3618, 2012.

[13] L. Goupil, E. Chanthery, L. Travé-Massuyès, and S. De-
liautier, “Tree based diagnosis enhanced with meta knowl-
edge,” in *34th International Workshop on Principles of
Diagnosis (DX’23)*, 2023.

[14] E. Real, C. Liang, D. So, and Q. Le, “AutoML-zero: Evolving machine learning algorithms from scratch,” in *International Conference on Machine Learning*, vol. 119. Pmlr, 13–18 Jul 2020, pp. 8007–8019.

[15] G. Pillonetto, F. Dinuzzo, T. Chen, G. De Nicolao, and L. Ljung, “Kernel methods in system identification, machine learning and function estimation: A survey,” *Automatica*, vol. 50, no. 3, pp. 657–682, 2014.

[16] J. Schoukens and L. Ljung, “Nonlinear system identification: A user-oriented road map,” *IEEE Control Systems Magazine*, vol. 39, no. 6, pp. 28–99, 2019.

[17] J. R. Koza, *On the programming of computers by means of natural selection*, ser. Koza, John R. Genetic programming. Cambridge, Mass. u.a.: MIT Press, 1992.

[18] M. Cranmer, “Interpretable machine learning for science with pysr and symbolic regression.jl,” 2023.

[19] M. S. Branicky, *Introduction to Hybrid Systems*. Boston, MA: Birkhäuser Boston, 2005, pp. 91–116.

[20] I. Saberi, F. Faghih, and F. S. Bavil, “A passive online technique for learning hybrid automata from input/output traces,” *ACM Transactions on Embedded Computing Systems*, vol. 22, no. 1, oct 2022.

[21] B. O. Bouamama, R. M. Alaoui, P. Taillibert, and M. Staroswiecki, “Diagnosis of a two-tank system,” In-tern Report of CHEM-project, USTL, Tech. Rep., 2001.

[22] N. Zaupa, L. Martínez-Salamero, C. Olalla, and L. Zaccarian, “Results on hybrid control of self-oscillating resonant converters,” in *IFAC Conference on Analysis and Design of Hybrid Systems (ADHS)*, vol. 54, 2021, pp. 211–216.

[23] ——, “Hybrid control of self-oscillating resonant converters,” *IEEE Transactions on Control Systems Technology*, vol. 31, no. 2, pp. 881–888, 2023.