📊 statistics

Identification of Gaussian Process State-Space Models with Particle Stochastic Approximation EM

This paper presents a maximum likelihood identification method for Gaussian process state-space models that combines a stochastic approximation EM algorithm with particle Markov chain Monte Carlo techniques to efficiently estimate parameters while preserving the full nonparametric description of system dynamics.

Original authors: Roger Frigola, Fredrik Lindsten, Thomas B. Schön, Carl E. Rasmussen

Published 2026-06-04

📖 5 min read🧠 Deep dive

CC BY 3.0

Original authors: Roger Frigola, Fredrik Lindsten, Thomas B. Schön, Carl E. Rasmussen

Original paper licensed under CC BY 3.0 (http://creativecommons.org/licenses/by/3.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to teach a robot how to walk, but you don't have a manual, and you don't know the laws of physics that govern its legs. You only have a video of the robot stumbling around. Your goal is to build a mathematical model that predicts where the robot will step next, but with a very special twist: the model must also tell you how unsure it is about that prediction.

This paper presents a new way to build that model using something called Gaussian Process State-Space Models (GP-SSMs). Here is a simple breakdown of what the authors did and why it matters.

The Problem: The "Black Box" of Motion

Usually, when scientists model a moving system (like a robot or a car), they try to guess the exact formula that describes how it moves. They might say, "It moves like a spring," or "It moves like a pendulum."

The Risk: If you guess the wrong formula, your model fails.
The Alternative: You could use a "generic" tool, like a neural network, to learn the shape of the movement. But these tools are like rigid molds; you have to decide exactly how many layers and knobs they have before you start. If you guess the size wrong, the model is either too simple to learn or too complicated and starts "hallucinating" (overfitting).

The Solution: The "Shape-Shifting" Clay

The authors propose using Gaussian Processes (GPs). Think of a GP not as a fixed formula, but as a lump of infinite, shape-shifting clay.

Instead of forcing the data into a pre-made mold, you let the clay naturally take the shape of the data.
The Magic: This clay has a special property: it knows how "smooth" or "jagged" the movement should be based on a few simple settings (called hyper-parameters).
The Result: You get a model that is flexible enough to capture complex, weird movements, but it doesn't need you to guess the exact structure of the math beforehand.

The Challenge: Finding the Settings

While the clay is great, you still need to tune the settings (the hyper-parameters) to make the model fit your data perfectly.

The Difficulty: The math to find the perfect settings is incredibly hard because the "clay" has infinite possibilities. It's like trying to find the perfect temperature for a cake when the oven has infinite temperature settings and you can't open the door to check.
The Paper's Trick: The authors use a method called Particle Stochastic Approximation EM (PSAEM).
- Imagine you are trying to find the best route through a foggy forest. You can't see the whole path.
- Instead of trying to map the whole forest at once, you send out a team of explorers (called particles).
- These explorers wander around, testing different paths.
- The algorithm listens to the explorers, learns from their mistakes, and slowly adjusts the map settings until the path becomes clear.
- Crucially, this method is very efficient. It doesn't need to check every single possibility; it just needs enough explorers to get a good guess.

What They Tested

The authors tested this method on two types of systems:

A Simple Linear System (The Straight Line):
- They fed the model data from a system that moved in a perfectly straight, predictable line.
- The Result: Even though the model was designed to handle complex curves, it successfully figured out that the movement was actually a straight line. It found the right settings and predicted the future steps accurately.
- The Confidence Check: When the model predicted a point far away from the data it had seen, it drew a wide error bar (a big "I'm not sure" zone). When it predicted a point close to the data, the error bar was narrow (a confident "I know this" zone). This is exactly what you want from a smart model.
A Complex Nonlinear System (The Wobbly Curve):
- They fed the model data from a system that moved in a very tricky, wobbly way with sharp turns.
- The Result: The model successfully learned the complex shape of the movement without being told what the shape was.
- The Confidence Check: In areas where the data was sparse or the movement was chaotic, the model correctly showed high uncertainty (wide error bars). In areas where it had seen similar movement before, it was confident.

The Big Takeaway

The main achievement of this paper is a method that allows computers to learn how systems move without needing a human to write down the physics equations first.

It's Flexible: It can learn simple straight lines or crazy, wobbly curves.
It's Honest: It doesn't just give you an answer; it tells you how confident it is in that answer.
It's Safe: By averaging over all possible shapes the movement could take, it avoids the trap of "overfitting" (memorizing the noise instead of learning the pattern).

In short, the authors built a tool that lets a computer "feel" its way through the physics of a system, learning the rules as it goes, while keeping a running tally of how sure it is about those rules. This is particularly useful for robotics and control systems where knowing what you don't know is just as important as knowing the answer.

Technical Summary: Identification of Gaussian Process State-Space Models with Particle Stochastic Approximation EM

Problem Statement
The paper addresses the challenge of identifying nonlinear dynamical systems where the underlying state transition dynamics are unknown and complex. Traditional system identification often relies on parametric models (e.g., specific equations of motion or fixed-structure neural networks). However, selecting an appropriate parametric structure is difficult without detailed prior knowledge, and fixing the number of parameters inherently bounds the complexity of the functions that can be fitted, potentially leading to underfitting or overfitting.

The authors propose using Gaussian Process State-Space Models (GP-SSMs) to overcome these limitations. A GP-SSM treats the state transition function as a random variable drawn from a Gaussian Process (GP) prior, providing a Bayesian nonparametric representation. This allows the model to capture a wide range of nonlinear dynamics while quantifying uncertainty. The specific problem tackled is the maximum likelihood (ML) identification of the hyper-parameters governing these GP priors and the noise characteristics, while retaining the full nonparametric description of the dynamics.

Methodology
The core difficulty in ML estimation for GP-SSMs is that the likelihood function requires marginalizing (integrating out) the latent state trajectories ( $x_{0:T}$ ) and the unknown transition functions ( $f_{1:T}$ ). While the unknown functions can be marginalized analytically to yield a prior over state trajectories, the resulting integral over the state space is analytically intractable.

To solve this, the authors employ the Particle Stochastic Approximation Expectation-Maximization (PSAEM) algorithm. This approach combines two key components:

Stochastic Approximation EM (SAEM): Instead of computing the exact expectation of the complete data log-likelihood (which is intractable), SAEM maintains a stochastic approximation of the auxiliary $Q$ -function. This is updated iteratively using samples from the smoothing distribution.
Particle Markov Chain Monte Carlo (PMCMC): Since direct sampling from the smoothing distribution $p(x_{0:T} | y_{0:T}, \theta)$ is impossible, the authors utilize Particle Gibbs with Ancestor Sampling (PGAS). PGAS constructs a Markov kernel that leaves the smoothing distribution invariant, allowing for the generation of sample trajectories $x_{0:T}^{[k]}$ to update the $Q$ -function approximation.

The algorithm proceeds iteratively:

E-step (Simulation): A new state trajectory is sampled using PGAS given the current parameter estimates.
Update: The stochastic approximation of the $Q$ -function is updated using the log-likelihood of the new trajectory.
M-step (Optimization): The parameters $\theta$ (hyper-parameters of the GP mean/covariance and noise variances) are updated by maximizing the approximated $Q$ -function using a quasi-Newton method (BFGS).

Key Contributions

Algorithmic Framework: The paper presents a method for maximum likelihood identification of GP-SSMs that avoids the need for pseudo-data sets or weak labels of unobserved states, which were limitations in prior approaches.
Integration of PGAS and SAEM: It demonstrates the application of Particle Gibbs with Ancestor Sampling within a Stochastic Approximation EM framework to handle the intractable integration required for GP-SSM identification.
Nonparametric Richness: The method retains the full nonparametric nature of the system model. By analytically marginalizing the nonparametric part of the model, the approach effectively averages over all possible dynamics consistent with the GP prior and the data, acting as a safeguard against overfitting (similar to Empirical Bayes).

Experimental Results
The authors validate the method on two types of systems:

Linear System: A linear dynamical system was identified using both linear and squared exponential covariance functions.
- When using a linear covariance function, the model successfully converged to parameters equivalent to a linear state-space model with Gaussian priors.
- When using a squared exponential (nonlinear) covariance function on linear data, the model correctly identified the system but exhibited high uncertainty (wide error bars) in regions where the data was sparse or the model was less confident, demonstrating the model's ability to quantify uncertainty even when the true system is simpler than the model class.
Nonlinear System: A complex nonlinear system with a sharp gradient in the state transition function and a quadratic measurement function (inducing bimodality) was identified.
- Using a Matérn covariance function (which allows for less smoothness than squared exponential) in the state direction, the model successfully captured the sharp transitions and the underlying dynamics.
- The identified model provided probabilistic predictions with error bars that reflected the uncertainty, particularly in regions far from the training data.

Significance and Claims
The paper claims that GP-SSMs offer a flexible framework for nonlinear system identification that incorporates high-level assumptions (like smoothness) without restricting the model to a specific parametric form. The proposed PSAEM approach enables the estimation of model hyper-parameters via maximum likelihood while preserving the nonparametric richness of the dynamics.

The authors emphasize that this approach provides probabilistic predictions with quantified uncertainty, which is crucial for applications like adaptive control and robotics where control strategies may depend on the confidence in the system dynamics. The method is presented as a robust alternative to fully Bayesian solutions (which require priors on parameters) and standard parametric identification, offering a balance between flexibility and computational feasibility through the use of particle-based inference.