Willems' Fundamental Lemma with Large Noisy Fragmented Dataset

Abstract: Willems' Fundamental Lemma enables parameterizing all trajectories generated by a Linear Time-Invariant (LTI) system directly from data. However, this lemma relies on the assumption of noiseless measurements. In this paper, we provide an approach that enables the applicability of Willems' Fundamental Lemma with a large noisy-input, noisy-output fragmented dataset, without requiring prior knowledge of the noise distribution. We introduce a computationally tractable and lightweight algorithm that, despite processing a large dataset, executes in the order of seconds to estimate the invariants of the underlying system, which is obscured by noise. The simulation results demonstrate the effectiveness of the proposed method.

• The paper introduces a robust algorithm that extends Willems’ Fundamental Lemma to large, noisy, and fragmented LTI datasets using ensemble averaging and invariant subspace recovery via SVD.
• It leverages statistical principles and grid search to accurately estimate unknown noise moments, achieving convergence at a rate of O(1/Nt) with precise rank diagnostics.
• Simulation results on synthetic systems validate the method’s high noise moment estimation accuracy, low computational cost, and practical scalability for data-driven control applications.

Extension of Willems' Fundamental Lemma to Large Noisy Fragmented Datasets

Motivation and Context

Willems' Fundamental Lemma (WFL) serves as a cornerstone for non-parametric, data-driven analysis and control of linear time-invariant (LTI) systems. The lemma allows any trajectory of the system to be described as a linear combination of data-driven trajectories, provided that the input data is persistently exciting. However, WFL assumes access to noiseless, contiguous trajectories—a highly restrictive condition in practical scenarios where input and output measurements are typically corrupted by noise, and where datasets are often collected over fragmented or interrupted experiments with missing intervals.

Traditional approaches to mitigating measurement noise in behavioral systems typically rely on low-rank regularization or robust control via uncertainty sets, leading to computational inefficiency or conservatism, especially in worst-case robust designs. The present work addresses these limitations by developing an algorithm for invariant subspace estimation tailored to datasets characterized by noise, fragmentation, and large scale, without any prior information on noise distribution or parameters. This approach capitalizes on the law of large numbers (LLN) across ensembles of experiments, enabling the reconstruction of system invariants under minimal a priori assumptions.

Theoretical Foundations

Multi-Experiment Setting

The formulation operates in a generalized multi-experiment context, where one collects a dataset $\mathcal{D}$ comprising $N_t$ disjoint experiments. Each experiment contains input-output sequences of fixed length, which may arise from splitting single trajectories or aggregating separate recording sessions. The only critical assumption is that collective input sequences across all fragments are persistently exciting of a sufficient order. The key insight is that system-invariant subspaces, specifically the left null spaces of the stacked Hankel matrices, remain shared across all fragments despite noise and fragmentation.

Noise Model and Empirical Estimation

Both input and output measurements are assumed contaminated by independent, zero-mean, additive i.i.d. noise, with all moments up to fourth order being finite but otherwise unknown. The theoretical aim is to estimate the left null space $\mathcal{N}(\mathbb{H}^{(i)})$ associated with noiseless Hankel matrices from only noisy, fragmented data.

Central to the method is constructing an empirical estimator matrix, $\hat{\mathbf{M}}_{\mathcal{D}}$, from ensemble-averaged empirical cross-correlations of noisy Hankel matrix rows, corrected by candidate noise moment values. The key result, established by Theorem 1, is that as $N_t \to \infty$, this empirical estimator converges (almost surely and at rate $\mathcal{O}(1/N_t)$) to the true aggregate correlation matrix (i.e., the expected best possible noise-free value), justifying estimation of invariant subspaces via singular value decomposition (SVD) of $\hat{\mathbf{M}}_{\mathcal{D}}$.

Algorithmic Solution

Parameter Estimation via Grid Search

Because the noise moments $m_1, m_2$ are strictly unknown, the approach involves a grid search across candidate moment parameters. For each candidate pair, the SVD of $\hat{\mathbf{M}}_{\mathcal{D}}(m_1, m_2)$ is computed, and the minimal singular value is tracked. The correct parameter pair is characterized by the lowest singular value, reflecting the cleanest recovery of the invariant left null space. This two-stage algorithm—empirical aggregation followed by grid-based parameter search and SVD—is designed to be computationally lightweight: its complexity is independent of dataset size, scaling only with the number of grid points and the (fixed) Hankel matrix block size.

Figure 1: Heatmap of the minimum singular value ($\sigma_{\min}$) over the search grid; optimal estimated noise moments sharply minimize $N_t$0.

Feasibility and Rank Diagnostics

In addition to examining singular values, the numeric rank of $N_t$1 is evaluated to confirm that the algorithm correctly recovers the theoretically anticipated deficiency (from Willems' Lemma: rank $N_t$2 in the noiseless case for $N_t$3 block rows), and solutions away from this basin are rejected.

Figure 2: Heatmap of the numerical rank of $N_t$4 across the same parameter search grid; the low-rank basin identifies correct invariant recovery.

Subspace Convergence and Statistical Guarantees

The geometric alignment of the recovered invariant subspace with the true null space (computed from exact, noiseless system matrices) is monitored as a function of the number of datasets $N_t$5, using principal subspace angles as the metric. The logarithmic decay of the maximum subspace error angle empirically validates the sharp consistency guarantee derived from the law of large numbers.

Figure 3: Logarithmic convergence of the maximum subspace error angle as sample size increases, demonstrating almost sure convergence of estimated system invariants.

Numerical Results

Simulations were conducted on a synthetic LTI system ($N_t$6), with $N_t$7 experiments, noisy input and output measurements, and a substantial noise bias ($N_t$8, $N_t$9). Across an exhaustive grid search ($\mathcal{N}(\mathbb{H}^{(i)})$0), the true noise moments were estimated with high accuracy ($\mathcal{N}(\mathbb{H}^{(i)})$1, $\mathcal{N}(\mathbb{H}^{(i)})$2). The algorithm’s empirical runtime was notably low (aggregation: 3.47s; search/SVD: 0.51s), confirming practical scalability. Subspace error decayed monotonically with increasing $\mathcal{N}(\mathbb{H}^{(i)})$3, substantiating robustness and theoretical guarantees.

Implications and Future Directions

This work establishes a non-parametric, robust approach for system identification and trajectory parameterization in realistic, large-scale, and fragmented experimental settings. By decoupling sample complexity from computational cost, the algorithm is compatible with emerging experimental and industrial scenarios characterized by parallelized data collection, incomplete datasets, and uncertain noise statistics. The general framework naturally accommodates extension to distinct input and output noise statistics (with grid search in higher dimensions), and provides a foundation for non-asymptotic statistical error bounds via concentration inequalities.

Practical implications include enhanced reliability of data-driven control design workflows in settings with fragmented measurements (e.g., networked control, distributed sensing, or experimental fault tolerance), with immediate relevance for robust data-enabled model predictive control, behavioral system identification, and learning-based stabilization. At a theoretical level, this work reinforces the connection between the law of large numbers and geometric system theory, strengthening the interface between modern data science and classical control.

Future directions include analysis of the trade-off between grid granularity, estimation bias for finite sample regimes, and extension to output feedback or input multiplicative noise models. Additionally, integrating Bayesian or variational approaches to estimate noise moments could further increase efficiency for high-dimensional systems.

Conclusion

The paper rigorously extends Willems' Fundamental Lemma to settings with large, noisy, and fragmented experimental data, providing provable, consistent, and efficient estimation of system invariants without requiring any prior knowledge of noise statistics. The methodology leverages ensemble averaging, SVD-based null space recovery, and effective computational strategies, with strong guarantees on both statistical consistency and numerical scalability.