Functional observability and subspace reconstruction in nonlinear systems

Time-series analysis is fundamental for modeling and predicting dynamical behaviors from time-ordered data, with applications in many disciplines such as physics, biology, finance, and engineering. Measured time-series data, however, are often low dimensional or even univariate, thus requiring embedding methods to reconstruct the original system's state space. The observability of a system establishes fundamental conditions under which such reconstruction is possible. However, complete observability is too restrictive in applications where reconstructing the entire state space is not necessary and only a specific subspace is relevant. Here, we establish the theoretic condition to reconstruct a nonlinear functional of state variables from measurement processes, generalizing the concept of functional observability to nonlinear systems. When the functional observability condition holds, we show how to construct a map from the embedding space to the desired functional of state variables, characterizing the quality of such reconstruction. The theoretical results are then illustrated numerically using chaotic systems with contrasting observability properties. By exploring the presence of functionally unobservable regions in embedded attractors, we also apply our theory for the early warning of seizure-like events in simulated and empirical data. The studies demonstrate that the proposed functional observability condition can be assessed a priori to guide time-series analysis and experimental design for the dynamical characterization of complex systems.


I. INTRODUCTION
Reconstructing the state space of complex dynamical systems is a key step for their quantitative understanding and forecasting. To do so, time-series analysis methods have been developed to characterize and model dynamical behaviors from recorded time-ordered data. In practice, such methods are constrained by the available measurement processes and data quality: data are often irregularly sampled, noisy, relatively short, and univariate. In cases that measured time series are expected to be lower dimensional compared to the system dynamics, the original state space can be reconstructed with embedding methods [1][2][3][4]. Embedding has been successfully applied across many fields, including for the characterization of chaotic dynamics [5,6], ecological and economic modeling [7,8], financial forecasting [9,10], medical diagnosis [11,12], and detection of dynamical transitions in palaeoclimate data [13] and oil-water flows [14]. These applications assume that an embedding is possible, that is, that a map from the time-series data to the original system's state space exists. However, such assumption may not hold in general: an established link between embedding and observability theories shows that reconstructing the entire original system is not always possible or suitable depending on the available time series [15,16].
The a priori conditions under which the entire system * arthur.montanari@uni.lu state can be reconstructed from available measurements are determined by the observability property [17]. In the linear case, observable systems satisfy the necessary and sufficient conditions for the inference of the full-state of the system, for example, via embedding methods and state estimators like Luenberger observers [18] or Kalman filters [19]. For nonlinear systems, a generalized notion of observability [20] sets a sufficient condition for the existence of an invertible mapping (diffeomorphism) between the system's original state space and the differential embedding space constructed from a given measurement function [15,16]. Not only the relation between embedding and observability determines the possibility of reconstructing the original dynamical system from the embedding of timeseries data, but it also characterizes how good such reconstruction is. For example, the embedding space of the Rössler system was empirically shown to have singularity points depending on the measured variable [21], which was later theoretically proven to be a consequence of unobservable regions in the original space [22]. Another case is the well-known Lorenz system: Fig. 1a-c shows that unobservable regions hamper the quality of the attractor reconstruction from an embedded time series. Consequently, the applicability and performance of methods based on embedding [23][24][25] and state estimation [26][27][28] are highly dependent on the observability properties of the system. Building from chaotic systems, often used as benchmark cases due to their short horizon of predictability, observability studies showed the potential to foster further discoveries in complex systems, (c) Reconstruction of the original coordinates of the Lorenz system from the differential embedding. The system is unobservable at the shaded plane (x1 = 0) and, therefore, the reconstructed attractor (left side) shows large errors closer to this region compared to the ground-truth (unmeasured) attractor (right side). (d) Reconstruction of the low-dimensional subspace (x1, x2) of the Lorenz attractor. The system is functionally observable everywhere with respect to the functional z = g(x) = x2 and, therefore, the reconstructed subspace (left side) is accurate compared to the ground-truth subspace (right side). Brighter colors correspond to states of the attractor closer to the system's unobservable plane.
For many systems and applications, though, complete observability is a condition that may be too restrictive. In practice, even if the original state space is not entirely observable (reconstructible), one may focus on particular subspaces (e.g., state variables) that are relevant to the considered applications. Examples include the estimation of the phase variable of nonlinear oscillators for synchronization analysis of chaotic systems [21,39,40], modeling of climate dynamics [41], and forecasting of financial crashes [42]; the positioning and tracking of a particular spatial coordinate (e.g., altitude) in autonomous aerial vehicles from indirect measurements [43]; or the inference of control variables (which dictate how close a system is to a bifurcation) for the early warning of transitions from healthy to disease states in atrial fibrillation [44] and epileptic seizures [45].
These practical problems motivate the concept of functional observability [46,47], which establishes conditions under which a desired functional of the system variables can be inferred from the available measurements (e.g., via the design of functional observers [46,48,49]). However, in spite of several applications designed for feedback control [49,50], fault detection [51], and, more recently, large-scale networks [52], the functional observability property is still restricted to linear dynamical systems. Therefore, a theory to establish conditions for the reconstruction of a nonlinear functional of a nonlinear system, and that provides guidance on how to perform and leverage this reconstruction, is still missing.
In this paper, we provide a generalization of the functional observability property to nonlinear dynamical sys-tems. This establishes a sufficient condition for the reconstruction of a nonlinear functional of state variables from a (possibly nonlinear) measurement function. If the system is functionally observable, we show how to determine the mapping between the differential embedding space and the original functional sought to be reconstructed, also proposing a coefficient of functional observability to locally characterize the quality of the functional reconstruction. Fig. 1 shows that, even though a system may not be completely observable (thus hampering the state-space reconstruction), it might still be functionally observable with respect to some low-dimensional subspace of interest, in which accurate reconstruction is still feasible. To illustrate the theoretical advantages of such framework in interpreting the effects of singularities and symmetries in the embedded state space, we present numerical simulations for chaotic benchmark systems with contrasting observability properties. Finally, we apply our theory for the analysis of a phenomenological model of seizure-like events, known as Epileptor [45]. We demonstrate that the presence of unobservable regions in the Epileptor's attractor can be used to provide early-warning signals of transitions from normal to seizure states in both simulated and empirical data.
The paper is organized as follows. Section II provides a background on the complete observability of dynamical systems. Section III presents the theoretical results for the generalization of functional observability to nonlinear systems. Section IV presents and discusses the numerical results in chaotic benchmark systems. Section V applies the proposed framework for the analysis of the Epileptor model and early warning of seizures. Finally, Section VI concludes the work.

II. BACKGROUND ON OBSERVABILITY OF NONLINEAR SYSTEMS
This section provides background on observability theory for nonlinear systems [20,38,53]. Consider the following nonlinear dynamical system where x ∈ X ⊆ R n is the state vector, y ∈ R q is the output vector (measurements), and f : X → V(X ) ⊆ R n and h : X → H(X ) ⊆ R q are smooth nonlinear functions. The explicit dependence of time in x(t) and y(t) is often omitted throughout this work. Let the flow map Φ T (x(t 0 )) : X → X be the solution of (1): The notion of local observability is formalized as follows.
Definition 1. [20,53] The nonlinear system (1), or the pair {f , h}, is locally observable at x 0 if there exists a neighborhood U ⊆ X of x 0 such that, for every state The system is said to be locally observable if it is locally observable at every x 0 ∈ X .
Definition 1 states that a system is locally observable around an initial state x(t 0 ) if x(t 0 ) can be uniquely reconstructed from the measurements y over a finite trajectory, as defined by the composition map h•Φ T (x(t 0 )). The local observability of a nonlinear system can be verified through the following algebraic condition.

Definition 2.
[53] The observable space O(x) of system (1) is the linear space of functions over the field R spanned by all functions of the form where L ν f h j (x) denotes the ν-th Lie derivative of the j-th component of h(x) along the vector field f (x), and s is the smallest integer such that ∇L k f h j (x) belongs to the span formed by functions (3) for all k > s. By definition, and where ∇ is the gradient operator with respect to x.
Theorem 1. The system (1), or the pair {f , h}, is locally observable at x 0 if there exists a neighborhood U ⊆ X of x 0 such that holds for every x ∈ U ⊆ X , where O(x) is the observable space.
Theorem 1 is a generalization of the observability property of linear systems to nonlinear systems. If the pair {f , h} is given by the linear functions Ax and Cx, where A ∈ R n×n and C ∈ R q×n , then it follows that the observable space O(x) is defined by the row space of Kalman's observability matrix [17,55] and therefore the linear system is observable if and only if rank(∇O) = n.

A. Observability and embedding
Differential embedding relates the original coordinates of a dynamical system to the derivatives of the measured time series. Formally, given some measurement y = h(x) over time t ∈ [0, T ], a differential embedding space E can be constructed via an appropriate choice of higher-order derivatives of y as coordinates: where y (ν) j denotes the ν-th time derivative of the jth component of the measured signal y. If the map Ψ : X → E is a diffeomorphism, then the original state space X and the embedding space E are related by a smooth and invertible change of coordinates. Therefore, in applications where the original state space is not available (due to unmeasured state variables), the underlying dynamical system can be assessed via an embedding of the measured time series. In practice, measured time-series data are available in discrete time and hence the embedding can be constructed by either computing their differential coordinates y (ν) j (e.g., via regularization methods to denoise the derivatives [56,57]) or using timedelay coordinates y j (t − kτ ), for some time delay τ and k ∈ N, following Takens' theorem [58,59].
From Definition 1, {f , h} is observable if the map is locally left invertible (injective) for some 0 ≤ ν ≤ s, That is, if x is uniquely determinable from y and its successive derivatives. The map (8) is equivalent to the set of functions (3) spanning the observable space O(x) [15], hence Ψ : X → E ≡ O(x). This equivalence establishes a direct relation between the theories of observability and embedding [15,16]: following the Inverse Function Theorem [53, Section 7.1], Ψ(x) is locally invertible if its Jacobian matrix has full (column) rank (i.e., if condition (5) holds). Therefore, a state x 0 is only reconstructable from the differential embedding of its measurements y over a finite time interval if Ψ is locally invertible or, in other words, the system is locally observable at x 0 .
Note that the choice of differential coordinates y (ν) j for the embedding space E is not unique and that an inappropriate selection of linearly dependent derivatives of y (i.e., functions of form (3)) may lead to dim(∇E) < dim(∇O(x)) [16]. Accordingly, we assume throughout this work that the embedding space (7) is defined by a minimum selection of linearly independent functions (3) such that dim(∇E) = dim(∇O(x)) around some neighborhood of x 0 . For example, we consider the embedding space E = {y,ẏ,ÿ} for the reconstruction of the Lorenz attractor ( Fig. 1a-c).

III. FUNCTIONAL OBSERVABILITY OF NONLINEAR SYSTEMS
Complete observability characterizes the sufficient condition for the (local) reconstruction of the full-state vector x of a dynamical system (1) from measurements y over finite time. Nonetheless, reconstructing the entire state vector x is often unfeasible or unnecessary in practice, and only a lower-dimensional function or subspace might be of interest, defined as where z ∈ R r is the vector sought to be reconstructed, often with dimension r n, and g : X → G(X ) ⊆ R r is a nonlinear smooth functional.
In what follows, we provide a generalization of the observability property [20], termed functional observability. Given a nonlinear dynamical system (1), functional observability establishes the conditions under which the functional (9) is reconstructible from the measured signal y, without necessarily requiring the full state x to be reconstructible. Therefore, a system may be functionally observable with respect to some functional (9) even though it is (completely) unobservable. Our results also generalize the functional observability property, originally established for linear systems [47], to dynamical systems described by a nonlinear vector field f , a nonlinear measurement function h and a nonlinear functional g. Fig. 2 summarizes the relation between our theory and previous works.
Before presenting the main result, we first formally define functional observability as a generalization of complete observability (Definition 1): Observability of linear systems (Kalman, 1959) Functional obsv. of linear systems (Jennings et al., 2011) Observability of nonlinear systems (Hermann and Krener, 1977) (1) and (9), or the triple {f , h, g}, is locally functionally observable at x 0 if there exists a neighborhood U ⊆ X of x 0 such that, for every state g( The system is said to be locally functionally observable if it is locally functionally observable at every x 0 ∈ X . Analogous to Definition 1, Definition 3 states that a system is locally functionally observable around an initial state x(t 0 ) if the functional g(x(t 0 )) can be uniquely reconstructed from the measurements y over a finite trajectory. Analogous to Definition 2, we define the functional space related to (9) as follows: (1) and (9) is the linear space of functions over the field R spanned by all functions of the form where µ is the smallest integer such that ∇L k f g j (x) belongs to the span formed by (4) for all k > µ.
Based on Definitions 2 and 4, we now establish the condition for the functional observability analysis of nonlinear systems. Consider a locally unobservable system at x 0 , i.e., Recall the theorem for the decomposition of unobservable systems [53, Theorem 97]: there exists a diffeomorphism T on U such that choosing an appropriate state transformationx = T (x) yields the partitioned vector wherex a ∈ R k andx b ∈ R n−k correspond, respectively, to the observable and unobservable parts of the system in some neighborhood of T (x 0 ). The transformed vector field is now given bỹ wheref a : R k → R k andf b : R n → R n−k , and the transformed measurement function is given bỹ which will only depend onx a . After this change of coordinates, the system (1) and (9) is given by We now note that if the following condition holds locally for every x 0 ∈ U ⊆ X , then, for somex = for all x ∈ U. Consequently, the system is functionally observable given thatg depends only onx a , which is the state vector corresponding to the observable subsystem. Condition (16) provides a sufficient condition for the local functional observability of the nonlinear system (1) and (9), or the triple {f , h, g}, at x 0 ∈ U. Note that this condition is locally equivalent to F(x 0 ) ⊆ O(x 0 ), that is, that the functional space (the subspace to be reconstructed) should be contained inside the observable space (the subspace reconstructible from the measurement function). The following theorem provides a condition for the functional observability of a nonlinear system that is equivalent to condition (16), but easier to implement since it does not require calculation of the functional space F(x 0 ) and the involved higher-order Lie derivatives of g.
Theorem 2. The nonlinear system (1) and (9), or the triple {f , h, g}, is locally functionally observable at The functional observability condition (18) establishes an easy-to-implement test that can be directly verified on the system's equations, that is, functions {f , h, g}. We note that, even though its derivation is based on the system decomposition (15) presented in [53,Theorem 97], this condition does not require any system transformation or prior knowledge of its equivalence transformation map T (x). This provides a significant advantage to analyze the limitations in subspace reconstruction of dynamical systems where deriving the transformation map T (x) and its inverse can be computationally intensive, such as the Epileptor model studied in Section V.
If the triple {f , h, g} is given by linear functions Ax, Cx and F x, where A ∈ R n×n , C ∈ R q×n and F ∈ R r×n , then the observable space O(x) is defined by the row space of the observability matrix (6) and the functional observability of a linear system can be verified using the following rank condition derived in [47,60]: If the full-state vector is sought to be reconstructed (i.e., g(x) = x), then dim{∇g(x)} = n and, therefore, the functional observability condition (18) reduces to the complete observability condition (5). Likewise, in the linear case, complete observability is a special case of functional observability by considering F = I n , where I n is the identity matrix of size n, which reduces condition (19) to the classical condition rank(∇O) = n.
As an illustrative example, consider a dynamical system (1) and (9), or equivalently the triple {f , h, g}, defined by is the region of analysis. Following Theorem 2, we first need to determine a basis for the observable and functional spaces according to Definitions 2 and 4: Rows {3, 4, 5, ...} of ∇O(x) are linear combinations of the first and second rows, leading to dim{∇O} = 2, ∀x ∈ X . This shows that, according to Theorem 1, the pair {f , h} is not (completely) observable, i.e., it is not possible to reconstruct the entire state vector x from y = x 2 . However, the triple {f , h, g} may still be functionally observable. By inspection, it is easy to see that ∇g(x) is a linear combination of the rows of ∇O(x), therefore satisfying condition (18) of Theorem 2 for all x ∈ X . This example illustrates that, even though a system is locally unobservable, it can still be locally functionally observable.

A. Functional observability and embedding
Complete observability establishes a sufficient condition for the existence of the (local) left-inverse map Ψ −1 : E → X from an embedding space to the original state space. Here, we generalize this relation by showing that functional observability establishes a sufficient condition for the existence of a map Φ : E → G(X ) from the embedding space to the subspace sought to be reconstructed. Furthermore, we demonstrate how to construct such a map if the system is functionally observable.
Let {f , h, g} be a functionally observable system with an observable space O(x) of local dimension (11). Following [53,Theorem 97], there exists a diffeomorphism T on U such that the transformationx = T (x) partitions the state vector as in (12). Consequently, the triple {f , h, g} can now be represented by {f ,h,g} as in (15). The diffeomorphism T is not unique and can be designed by partitioning it as for some neighborhood U ⊆ X of x 0 . The functions Ψ a (x) and Ψ b (x) can be designed as follows: 1. Construct a map Ψ a (x) by selecting a minimum set of linearly independent functions of form (3) 2. Construct an arbitrary map Ψ b (x) such that dim{∇T (x)} = n, ∀x ∈ U ⊆ X , i.e., ∇T (x) has full rank.
Remind that Ψ a (x) and Ψ b (x) are only valid around a local neighborhood U ⊆ X of some state x 0 ∈ X . The transformation Ψ a (x) : U → E defines a basis for the embedding space E, which depends only on y and its successive derivatives (see Eq. (7)). Since y = h(x) = h(x a ) is a function only ofx a (following Eq. (14)), then step 1 guarantees thatx a ∈ E.
Step 2 constructs an arbitrary function Ψ b (x) : U → E c that defines a basis to the complement of the embedding space E c in order to accomplish the Inverse Function Theorem. Therefore, the designed diffeomorphism T : U → E ∪ E c has a local inverse map T −1 .
Since the system is functionally observable, then relation (17) holds and z = g(x) =g(x a ), which depends only onx a ∈ E. Therefore, there exists a map g : E →G(E) from the embedding space E to the functional sought to be reconstructed z =g(x a ). Equivalently, z = g(x) can be reconstructed from the composition map Φ = g • T −1 , which can be locally constructed from the known function g and the designed transformation T according steps 1 and 2. Finally, if the system is functionally observable, then   3 illustrates a commutative diagram of the composition map (24). Clearly, if k = n, then the relation between complete observability and embedding follows as a special case: Building on the previous example, given that the system (20) is functionally observable, we demonstrate how to compute (reconstruct) the sought vector z from the measurement signal y. To this end, we first decompose the system as in (12) by constructing a diffeomorphism T (x) partitioned as (23): wherex a ∈ R 2 ,x b ∈ R 1 , and dim{∇O} = k = 2. Note that Ψ a (x) = [yẏ] T defines a map between the observable vector and the differential embedding coordinates, while Ψ b (x) is chosen arbitrarily to accomplish the Inverse Function Theorem. Therefore, the diffeomorphism has the inverse function Consequently, z can be computed as As expected, we have thatg(x a ,x b ) ≡g(x a ) and, therefore, g T −1 (x) depends only onx a , which is a function of y and its subsequent derivatives.

IV. OBSERVABILITY OF CHAOTIC SYSTEMS
We explore the functional observability property in different types of chaotic dynamical systems with contrasting observability properties, considering different measurement functions as well as functionals sought to be reconstructed. Following the theoretical conditions established in Section III, functional (or full-state) reconstruction is possible when the system is functionally (or completely) observable. Beyond this binary characterization of the system observability (i.e., either the system is or is not observable), we show that the performance of the reconstructed functional (full-state) vector is dependent on the proximity of the system state to functionally (completely) unobservable regions in the state space. The reconstruction errors are related to the sensitivity of the maps Ψ −1 and Φ to small perturbations (e.g., noise in the measured signals y(t)), and can be quantified by the absolute condition number of the inverse maps between the embedding coordinates and the reconstructed state: We address the condition numbers κ(Ψ −1 ) and κ(Φ) as the coefficients of complete and functional observability, respectively; a nomenclature that was previously adopted for κ(Ψ −1 ) in studies restricted to complete observability [16,22,35,38]. Results show that these coefficients can be employed to assess the quality of the (functional) state reconstruction as x(t) approaches unobservable regions: the larger κ, the higher the reconstruction error in the corresponding states.
In what follows, chaotic systems were numerically integrated using a fourth-order Runge-Kutta integrator with time step dt = 0.01s for a total simulation time T = 1100s, where the initial transient T trans = 1000s was discarded and initial conditions were randomly drawn from a normal distribution (i.e., x i (0) ∼ N (0, 1), i = 1, . . . , n). Codes are publicly available at https://github.com/montanariarthur/ NonlinearObservability. The symbolic construction of Lie derivatives (3) spanning the observable space O(x), as well as maps Ψ and Φ, is illustrated in these codes. Note that Ψ is composed by the minimum set of linearly independent functions (3) spanning O(x). Therefore, as a special case for q = 1 and ν = n − 1, Ψ(x) = O(x).

A. Lorenz system: observability and symmetry
The well-known Lorenz'63 system is given by where (R, σ, b) = (28, 10, 8/3) is a set of parameters that leads to a chaotic attractor (Fig. 1c, right). Here, we consider that only the state variable x 1 is available for measurement (i.e., y = h(x) = x 1 ). The observability of the pair {f , h} can thus be verified through the observability matrix Since det(∇O(x)) = −σ 2 x 1 , then, following condition (5), the system is locally observable at every state x ∈ R 3 except at Given the differential embedding coordinates E = {y,ẏ,ÿ} (Fig. 1b), the entire state vector x can be reconstructed by computing the inverse map O −1 (y) : E → X (where, as a special case, Ψ(x) = O(x)). Theoretically, existence of this map is not guaranteed only when the system is locally unobservable. In this example, the unobservable subspace corresponds to the exact region in the state space where x 1 = 0, which has null Lesbegue dimension. However, in practice, this map degenerates as the system trajectory approaches the neighborhood of x 1 = 0, corresponding to a gradual loss of system observability [22]. Fig. 4a illustrates the reconstruction error e x (t) = x(t) −x(t) , wherex = O −1 (y,ẏ,ÿ) is the reconstructed (estimated) state vector, considering a small additive noise to the measured time series: . Note that noise is largely amplified and the reconstruction error increases significantly as the system state approaches the unobservable region (x 1 (t) → 0).
Despite the unobservability at x 1 = 0, the system {f , h, g 1 } is always functionally observable with respect to the functional g 1 (x) = x 2 , where row(∇g(x)) ⊆ row(∇O(x)), ∀x ∈ R 3 . In practice, z = g 1 (x) can be reconstructed from the composition map (24) given bỹ Fig. 4b shows the reconstruction error e z (t) = z(t) −ẑ(t) , whereẑ = Φ(y,ẏ). As expected, since the system is functionally observable for all x ∈ R 3 , the reconstruction error e z (t) of the functional is not affected by the unobservable region x 1 = 0 and remains bounded (see also Fig. 1d). Fig. 5 presents the coefficients of observability for the Lorenz system. The coefficient of complete observability increases as x 1 (t) → 0, indicating a substantial increase of sensitivity of the local reconstruction map Ψ −1 (x) to small perturbations, as observed in the large reconstruction errors e x (t) for x 1 (t) → 0 in Fig. 4a. On the other hand, the coefficient of functional observability remains well-conditioned and constant throughout the entire attractor, which is supported by the insensitivity to noise in the reconstruction error e z (t) in Fig. 4b. These results demonstrate that these coefficients can provide proxy indicators of the "practical" consequences of the lack of (functional) observability of these systems as the state approaches (functionally) unobservable regions, where local reconstruction maps become highly sensitivity to noise and, therefore, fail to provide an accurate reconstruction of the original state space.
The Lorenz system is marked by a clear relation between observability and symmetry. Since h(x) = x 1 is directly measured and g(x) = x 2 is functionally observable, one can observe (reconstruct) the dynamics in the (x 1 , x 2 ) plane of the Lorenz attractor. The functional observability of this plane is directly related to the global invariance of the Lorenz attractor under the [22]. Given that x 3 is invariant under this symmetry, one can only distinguish which "wing" of the chaotic attractor the system state belongs to at a given time instant t by accurately observing variables x 1 and x 2 (Fig 1d). Therefore, the functional observability of the triple {f , h, g 1 } provides the necessary and sufficient information for this characterization based on the measured time series y(t). Moreover, it is evident that the lack of complete observability in the system {f , h, g 1 } is due to variable x 3 , which can be rigorously verified by noting that the functional observability condition (18) is only satisfied for a triple {f , h, g 2 }, g 2 (x) = x 3 , if x 1 = 0.

B. Cord system: fast and slow dynamics
The Cord system, a variation of the Lorenz'84 system, is given by [61]     ẋ where (a, b, F, G) = (0.258, 4.033, 8,1). The chaotic attractor is illustrated in Fig. 6a. The system dynamics is marked by two oscillation modes with a clear timescale separation [62]. Oscillations in the slow timescale can be approximately monitored by the "slow phase" variable θ s = x 1 (Fig. 6b), in which a full revolution of the system is completed every time the trajectory approximates the cord filament close to the origin (defining the Poincaré section P = {x : x 1 = 0,ẋ 1 > 0}) [40]. Oscillations in the fast timescale, on the other hand, can be monitored by the "fast phase" variable θ f = tan −1 (x 2 /x 3 ) (Fig. 6c).
Here, we consider the measured time series y = h(x) = x 2 and that the slow and fast phase variables are the functionals sought to be reconstructed, i.e., g 1 (x) = θ s and g 2 (x) = θ f . Fig. 7 shows the coefficients of observability for the Cord system. Full-state reconstruction of the Cord system from the measured time series is not possible for a considerable range of states in the system trajectory (Fig. 6a), defined by the plane and is expected to be ill-conditioned when the system state is close to the vicinity of this plane. The unobservable plane can be visualized in the (x 2 , x 3 ) section of the attractor in Fig. 7a.
Similarly to the full-state reconstruction problem, the reconstruction of the system's slow timescale (i.e., g 1 (x)) from time-series data of a state variable dominated by a fast timescale (e.g., x 2 ) is hampered by the lack of observability in a large subspace of the state space, as indicated by the regions with very large coefficients of functional observability in Fig. 7b. Contrariwise, reconstruction of the fast timescale (i.e., g 2 (x)) from x 2 is wellconditioned throughout the entire system trajectory, except as (x 2 , x 3 ) → (0, 0) (Fig. 7c, red circle). The lack of functional observability at this singularity region in the attractor is not surprising: it corresponds exactly to the region in which the fast phase variable θ f is a (locally) ill-defined function, and the fast phase "collapses", undergoing an inversion of its rotational direction [61].
This example illustrates that, even though the system may have a (relatively) large unobservable region X u ⊂ X , one may find that, even in this unobservable region, the system can still be functionally observable with respect to some functional g(x) aside from a significantly smaller subregion X fu ⊂ X u . In this example, the region of interest X is the Cord attractor A, the "completely" unobservable region X u is the 2-dimensional plane defined by (35), and the functionally unobservable region is the 1-dimensional line X fu = {(x 1 , 0, 0) | x ∈ A}.
C. Hindmarsh-Rose system: a neuron model Building up from the chaotic benchmarks, we now consider a phenomenological model of neuron dynamics given by the Hindmarsh-Rose (HR) model [63]: where x 1 is the membrane potential, x 2 is the fast recovery current, and x 3 is the slow adaptation current. Providing both a simplification of the biophysical Hodgkin-Huxley neuronal model and a generalization of the FitzHugh-Nagumo model, the HR model can reproduce a wide range of dynamical behaviors, including quiescence and (irregular) spiking and bursting [64]. Moreover, depending on the bifurcation parameters, this system can also shift to chaotic regimes, as investigated both computationally [64] and experimentally [65]. Here, we consider the set of parameters lying in the chaotic regime: (a, b, c, d, I, r, s, x R ) = (1, 3, 1, 5, 3.25, 0.001, 4, −8/5).
One might wonder if, despite the lack of complete observability at x 1 = 0, the system {f , h 2 , g i } is still functionally observable with respect to, for example, g 1 (x) = x 1 or g 2 (x) = x 3 . However, unlike the previous examples, the HR model remains locally unobservable at x 1 = 0 with respect to both functionals, as observed in the coefficients of functional observability shown in Fig. 8b,c. Nevertheless, note that the neighborhood of x 1 = 0 where the reconstruction map is (locally) ill-conditioned is substantially smaller for {f , h 2 , g 1 } compared compared to {f , h 2 , g 2 }. These results suggest that reconstruction of g 1 (x) is more reliable than g 2 (x) in the presence of small perturbations as x 1 → 0.  Examining the local maps yield The presence of the terms x 1 and x 2 1 in the denominator of Eqs. (38) and (39) elucidate the results shown in Fig. 8b,c. The sensitivity to small perturbations in the reconstruction of functional g 1 (x) is only inversely proportional to the distance between x 1 and the unobservable region, whereas the sensitivity of the reconstruction of g 2 (x) is inversely proportional to the quadratic of this distance-leading to a highly ill-conditioned map for |x 1 | 1. This theoretical (local) analysis is also supported by computing the reconstruction maps Φ : E → G(X ) and evaluating the corresponding reconstruction performance for each functional. Fig. 9 shows that, in the presence of small measurement noise v(t) ∼ N (0, 0.01), reconstruction of z 2 = g 2 (x) yields very poor results, with a high rootmean-square error (RMSE) of 0.4006, compared to the RMSE of 0.0242 for the reconstructed vectorẑ 1 = g 1 (x).
As in the Cord example, reconstruction of the slow timescale dynamics (g 2 (x) = x 3 in the HR model) from time-series data corresponding to a variable dominated by the fast timescale (h 2 (x) = x 2 ) is marked by the presence of unobservable regions which significantly hamper the quality of the reconstruction in the vicinity of these regions. On the other hand, measuring a variable dominated by the fast timescale can still provide accurate reconstruction of other variables dominated by the same timescale (g 1 (x) = x 1 ). This relation between the timescale separation and functional observability of a system, with respect to variables belonging to the same or different timescales than the measured variable, can be observed both in the Cord and HR models.

V. EARLY WARNING OF SEIZURES
Characterizing and predicting epileptic seizures are long-standing challenges in clinical neuroscience [66,67]. Accurate and interpretable methods for the prediction of seizure events will drastically improve epilepsy management, providing early warnings to alert patients or trigger interventions [68]. On top of black-box and datagreedy deep learning algorithms, dynamical-based topological analysis can concur in discovering universal routes to epilepsy and foster new methods for early warning, many of which can be based on the embedding of timeseries data [45,69]. We investigate the observability and embedding properties of such application by considering a dynamical model describing seizure dynamics in the brain. The model, termed Epileptor [45], involves bifurcation dynamics to reproduce resting, spiking, and bursting behaviors observed in electroencephalogram (EEG) signals, modeling the multiple timescale oscillations recorded during epileptic seizures. The Epilep- tor is defined by [45] where (r 1 , r 2 , I 1 , I 2 , γ) = (−1.6, 1, 3.1, 0.42, 0.01) are the system parameters, (τ 0 , τ 2 ) = (2857, 10) are the timescale constants, and the coupling functions are given by This model consists of three subsystems with different timescales: (x 1 , x 2 ) governs the system's oscillatory behavior, (x 4 , x 5 ) introduces the spikes and wave components typical in seizure-like events, and x 3 represents a slow permittivity variable that determines how close the system is to the seizure threshold. Due to the slow-fast timescale separation induced by τ 0 , x 3 is usually interpreted as a quasi-steady state parameter [45], enabling bifurcation analysis.

A. Functional observability analysis
Monitoring the permittivity variable x 3 provides an early-warning signal of a dynamical transition from normal to seizure states in the Epileptor model. Despite the phenomenological nature of the model, this permittivity variable is most likely related to slowly changing biophysical parameters (e.g., extracellular processes or ionic concentrations) [45]. Given that such parameters are hardly measurable in biomedical setups, we investigate whether it is possible to infer the permittivity variable (i.e., the functional z = g(x) = x 3 ) from more easily accessible measurements, such as EEG recordings of seizure-like events (modeled as the measurement signal y = h(x) = x 1 + x 4 due to its close resemblance to actual EEG data [45]). Fig. 10a illustrates the dynamics of the functional z(t) and output y(t). Under the assumption that the Epileptor is a proper representation of the underlying process, a functional observability analysis of model (40) can establish if it is feasible to reconstruct this functional and, therefore, provide an early-warning signal of seizure events from EEG data. Due to the model complexity, an analytical derivation of the functionally unobservable regions of the Epileptor model is hardly tractable. Instead, Fig. 10b,c presents the coefficients of functional observability of the triple {f , h, g} computed over the system's attractor. The system alternates between two regions of the attractor, the normal state and the seizure state, as x 3 crosses predetermined thresholds marking bifurcation points (i.e., points in the parameter space where qualitative changes in system dynamics occur). While the coefficients κ are fairly well-conditioned in the seizure region of the attractor (κ < 10 2 ), the normal region has relatively larger coefficients (κ ≈ 10 4 ) with two remarkable ill-conditioned singularities (κ > 10 7 ) highlighted in Fig. 10b. This indicates the presence of two functionally unobservable states in the normal region of the Epileptor's state space, one of them located exactly at the saddle-node bifurcation point from normal to seizure regime (x 3 ≈ 2.9,ẋ 3 < 0 [70]). Introducing linear additive process noise to the model (40) promotes a larger exploration of system's state space, uncovering other functionally unobservable singularities in the normal region (Fig. 10c), including a few in the seizure region. Nonetheless, the analysis remains qualitatively similar between the deterministic and stochastic systems: both show considerably larger values of κ in the normal region compared to the seizure region (see also Fig. 12b). Consequently, high errors are expected in the reconstruction of the permittivity variable from the measured signal y(t) during the normal regime of the Epileptor.

B. Early-warning signals and observability
At first, large coefficients of functional observability in the Epileptor's normal region indicate that reconstructing the slow permittivity variable from EEG data is particularly challenging. However, our analysis established an interesting relation between the Epileptor's observability and topological features: normal (seizure) regions of the attractor correspond to regions with large (small) coefficients of functional observability. This relation can be explored to develop early-warning indicators of seizure-like events in simulated and empirical data.
Typical early-warning signals of critical transitions studied in the literature, such as variance and autocorrelation, are computed from time-series data [71]. Evaluating the system's observability, on the other hand, requires prior knowledge of the system's equations (1), often absent for real-world systems. Theory states that unobservable regions in the state space are associated to the loss of dimension of the observable space (i.e., condition (18) does not hold). As a consequence, closer to unobservable regions, embedded trajectories squeeze into a small low-dimensional neighborhood due to the loss of diffeomorphism between the embedding space and the original state space [15]. This phenomenon is illustrated in Fig. 11 for the embedded attractors of different dynamical systems with poorly observable regions as well as real-world data. Such local topological feature can be assessed by monitoring the smallest singular value σ de computed from an embedded time series with embedding dimension d e (see Appendix B for details). As σ de → 0, the effective dimension of the embedded time series drops, implying that the diffeomorphism between the embedded and original attractors is not locally preserved (and, therefore, the system is locally unobservable). In what follows, we apply the coefficient σ de , hereby referred to as "time series-based singular value decomposition" (tSVD), as a proxy measure of the system's observability computed from time-series data, which, as we show next, has a high correlation with the coefficients of observability (Fig. 12d). increases (decreases) during normal (seizure) regimes of the Epileptor, the tSVD decreases (increases). The anticorrelation is confirmed by a Pearson's correlation index of ρ = −0.96 between both coefficients in logarithmic scale (Fig. 12d, top). As expected, when the Epileptor switches to the normal region, which is poorly functionally observable (κ > 10 4 ), the smallest singular value σ de tends to zero (σ de ≈ 10 −16 ), implying an effective loss of dimension of the embedded time series. For the stochastic model, the broader state-space exploration of this system yields a higher variation of κ and σ de . Nevertheless, the same anti-correlation between κ and σ de can be observed by smoothing the coefficients over a moving average window. On average, κ (σ de ) increases (decreases) during the normal regime of the Epileptor, yielding a Pearson's correlation index of ρ = −0.82 (Fig. 12d,  bottom).
The behavior of the tSVD measure is remarkably aligned with that of typical early-warning signals used to detect critical transitions from time-series data [73], including those for seizure warning from EEG data [74]. Indeed, there is a sharp increase of the tSVD close to the dynamical transition from normal to seizure state (Fig. 12c). This may be attributed to the Epileptor's unobservability at the saddle-node bifurcation point, as pinpointed in Section V A. Fig. 11c shows that the embedded trajectories squeeze to singularity point at the unobservable (bifurcation) point-a feature that can be further explored for early-warning detection of seizure events. This topological characteristic of the embedded attractor is present not only in the Epileptor model but also in real data, as shown in Fig. 11d for the embedding space constructed from human intracranial EEG data (public data available at [72]; sampling protocols and preliminary analysis are described in [67]). Further computing the tSVD in human EEG data provides interesting results as illustrated for a representative patient in Fig. 13: the seizure onset is often preceded by a decrease of σ de followed by a sharp increase close to the critical transition, a characteristic that may be explored for real-time monitoring and detection of seizure events. The same pattern was observed for all patients in the considered dataset, although a thorough statistical investigation of the tSVD as a early-warning signal of seizure events (and other critical transitions in complex systems) is left for future work.

VI. DISCUSSION
The established relation between observability and embedding theories opens a new research direction of special interest to (nonlinear) time-series analysis. Our theory formally determines the conditions for reconstructing the system state from time-series data, often lowdimensional or univariate. In fact, measuring every relevant state variable is in practice constrained by physical limitations or operational costs. Hence, indirect estimation of unmeasured variables is required for the observation of physical, biological, ecological, and other complex dynamical systems.
For applications which require reconstructing only a few key state variables or lower-dimensional subspaces, we formalize the notion of functional observability for nonlinear systems. Our results can provide a priori knowledge of the reconstruction limitations and embedding features, depending on the available time-series data. We show that, even if a system is not completely observable (reconstructible), it may still be functionally observable with respect to the variables or subspace of interest. This provides useful insights about the dynamical system's properties and can be used to guide experimental design and data-processing methods according to the investigated hypotheses and available measurement processes.
In the context of systems biology, observing system dynamics is often hampered by technical limitations that prevent the simultaneous measurement of multiple biophysical variables (e.g., multiple ion channels in single neurons). By identifying conditions for accurate inference of variables from time-series data, the presented functional observability analysis can thus guide experimental design. Consider, for example, the HR neuron model investigated in Section IV C. The fast recovery current and the slow adaptation current represented by the system's state variables are related to transport rates of fast (e.g., sodium and potassium) and slow (e.g., calcium) ion channels, respectively [65]. Our analysis reveals structural limitations in the HR neuron model that prevent an accurate inference of calcium flux from measures of sodium/potassium fluxes. The opposite, instead, seems feasible, given that the system is completely observable everywhere when inferring sodium/potassium flux from measures of calcium flux.
Likewise, in a biomedical context, evaluating the functional observability of the Epileptor model shows that reconstructing the slow permittivity variable from EEG data (y = x 1 + x 4 ) is complicated by the system's poor observability in the attractor's normal state (Section V). Contrariwise, independently measuring the state variables x 1 and x 4 (i.e., y = [x 1 x 4 ] T ) yields well-conditioned coefficients of observability throughout the entire attractor. This suggests that applying data pre-processing methods in EEG time series to uncouple the oscillatory behavior (modeled by x 1 ) from the spikes and wave components (modeled by x 4 ) may lead to better performance in reconstructing the permittivity variable for early-warning of seizures.
In addition to applications, the proposed theory opens new theoretical research directions for many disciplines. First, the Cord and HR neuron examples show interesting links between the functional observability of a system and its intrinsic timescales. In both cases, high reconstruction errors stem from estimating slow variables from measures of fast variables. Future works can formally explore this interesting relation, complementing the analysis for linear systems [75], by extending the notion of functional observability to (nonlinear) differential-algebraic systems of form where a strong timescale separation arises from a quasisteady-state assumption (ẋ 2 ≈ 0). Second, our analysis of the Epileptor model shows a potential relation between the system's observability and its bifurcation points. Whether the loss of observability close to critical transitions is a universal behavior or a particularity of the Epileptor remains to be investigated. In ecological networks, time-series data of variables that make the system completely observable often lead to earlier warning of critical transitions [32]. Our time-series-based coefficient tSVD, aside from indirectly quantifying observability, may also capture features related to the Central Limit Theorem [76] (that, close to bifurcation points, dynamical systems can be locally reduced to low-dimensional normal forms). However, it is still to be investigated whether our framework only applies to transitions induced by local bifurcations, or it can be extended to other types like boundary crisis involving chaotic attractors [77,78]. The potential use of tSVD for early-warning detection of critical transitions in complex systems, similar to other signals like increasing variance and autocorrelation [73,79], may lead to promising theoretical developments and applications. Third, our theory fosters data-driven methods for the automated construction of the embedding space, and its map Φ to the original system's attractor, in applications where analytical analysis of the model is untractable (e.g., due to unknown parameters or high dimensionality). This would thus extend previous works on automated embedding construction [80] and full system identification from embedding coordinates [81]. Although our application examples focus on univariate measurements and functionals, the theory is formalized for multivariate cases (q, r ≥ 1) and can be directly applied to determine the existence and conditioning of such map, assessing how good the reconstruction is expected to be (locally).
Finally, for the study of high-dimensional problems, our results call for extensions based on graph-theoretical conditions [31,37,52] or network motifs [35]. In fact, as the computation of Lie derivatives is particularly demanding for high-dimensional systems, scalable strategies have yet to be developed to investigate the functional observability of large-scale nonlinear networks.