anomaly-detection-dynamical-systems

name: anomaly-detection-dynamical-systems description: "Log-likelihood ratio based real-time anomaly detection for dynamical systems with theoretical characterization of error rates. Identifies system anomalies from multiple possible plant models for linear Gaussian systems. Use for: anomaly detection, fault diagnosis, dynamical system monitoring, log-likelihood ratio testing, industrial control applications. Activation: anomaly detection, log-likelihood ratio, dynamical systems, industrial control, fault detection."

Anomaly Detection for Dynamical Systems

Log-likelihood ratio based real-time anomaly detection with theoretical error rate characterization for linear Gaussian systems.

Overview

This methodology provides:

• Real-time anomaly detection for industrial control applications
• Theoretical characterization of detection error rates
• Multiple plant model framework for different anomaly types
• Observer design leveraging error rate characterization

Core Concepts

Anomaly Detection Framework

System operates in modes:
├── Nominal mode (no anomaly)
├── Anomaly type 1
├── Anomaly type 2
└── ...

Each mode → Different plant model → Different likelihood

Log-Likelihood Ratio Test

For each candidate model M_i:
LLR_i = log(p(y| M_i) / p(y| M_0))

Where:
- y: observations
- M_i: model with anomaly i
- M_0: nominal model

Decision: Select model with maximum LLR

Error Rate Characterization

False alarm rate: P(choose M_i | M_0 true)
Miss detection rate: P(choose M_0 | M_i true)

These are characterized theoretically for linear Gaussian systems.

Mathematical Framework

Linear Gaussian System

x_{t+1} = A_i x_t + B_i u_t + w_t,  w_t ~ N(0, Q_i)
y_t = C_i x_t + D_i u_t + v_t,      v_t ~ N(0, R_i)

Where i indicates the operating mode (nominal or anomaly)

Log-Likelihood Computation

For Kalman filter estimates (x̂_t, P_t):

Innovation: ν_t = y_t - C_i x̂_t - D_i u_t
Innovation covariance: S_t = C_i P_t C_i^T + R_i

Log-likelihood increment:
ΔLLR = -0.5 * [ν_t^T S_t^{-1} ν_t + log|S_t| + m*log(2π)]

Where m is measurement dimension

Error Rate Analysis

For two models M_i and M_j:

KL divergence: D_{ij} = E_{M_i}[log(p(y|M_i)/p(y|M_j))]

This determines error rates:
- Lower KL divergence → Higher confusion probability
- Higher KL divergence → Better separability

Algorithm

Step 1: Model Bank Setup

class ModelBank:
    def __init__(self):
        self.models = {
            'nominal': SystemModel(A0, B0, C0, D0, Q0, R0),
            'anomaly_1': SystemModel(A1, B1, C1, D1, Q1, R1),
            'anomaly_2': SystemModel(A2, B2, C2, D2, Q2, R2),
            # ...
        }
        self.kalman_filters = {
            name: KalmanFilter(model) 
            for name, model in self.models.items()
        }

Step 2: Real-time Detection

def detect_anomaly(observation, control_input):
    """
    Real-time anomaly detection using LLR.
    """
    llrs = {}
    
    for model_name, kf in kalman_filters.items():
        # Predict
        kf.predict(control_input)
        
        # Update with observation
        innovation = kf.compute_innovation(observation)
        
        # Compute log-likelihood increment
        ll_increment = compute_log_likelihood(innovation, kf.innovation_cov)
        
        # Accumulate LLR
        if model_name not in llrs:
            llrs[model_name] = 0
        llrs[model_name] += ll_increment
    
    # Select most likely model
    detected_model = max(llrs, key=llrs.get)
    confidence = compute_confidence(llrs)
    
    return detected_model, confidence, llrs

Step 3: Error Rate Estimation

def estimate_error_rates(model_i, model_j):
    """
    Theoretical error rate between two models.
    """
    # Compute steady-state Kalman gains
    K_i = steady_state_kalman_gain(model_i)
    K_j = steady_state_kalman_gain(model_j)
    
    # Compute asymptotic covariance
    Σ_i = asymptotic_innovation_covariance(model_i)
    Σ_j = asymptotic_innovation_covariance(model_j)
    
    # KL divergence approximation
    D_ij = 0.5 * (trace(Σ_j^{-1} Σ_i) - dim + log(det(Σ_j)/det(Σ_i)))
    
    # Error rate (Chernoff bound approximation)
    error_rate = exp(-D_ij / 2)
    
    return error_rate

Observer Design Leveraging Characterization

Design Criterion

Given: Required error rates (P_FA, P_MD)
Find: Observer parameters (L, Q_observer, R_observer)

Such that:
- P(false alarm) ≤ P_FA
- P(miss detection) ≤ P_MD

Optimization Problem

def design_observer(models, requirements):
    """
    Design observer to meet error rate requirements.
    """
    def objective(observer_params):
        observer = build_observer(observer_params)
        
        # Simulate error rates
        fa_rate = simulate_false_alarms(observer, models)
        md_rate = simulate_miss_detections(observer, models)
        
        # Penalty for constraint violation
        penalty = max(0, fa_rate - requirements['P_FA']) + \
                  max(0, md_rate - requirements['P_MD'])
        
        return penalty
    
    # Optimize observer parameters
    optimal_params = minimize(objective, initial_guess)
    return build_observer(optimal_params)

Use Cases

1. Industrial Process Monitoring

Application: Chemical plant fault detection
Models:
- Nominal: Normal operating conditions
- Anomaly 1: Sensor drift
- Anomaly 2: Valve malfunction
- Anomaly 3: Heat exchanger fouling

Benefit: Early detection enables prompt maintenance

2. Aerospace Systems

Application: Aircraft subsystem monitoring
Models:
- Nominal: Standard flight conditions
- Anomaly 1: Actuator degradation
- Anomaly 2: Sensor bias
- Anomaly 3: Structural damage indicators

Benefit: Safety-critical fault detection

3. Power Grid Monitoring

Application: Electrical grid anomaly detection
Models:
- Nominal: Stable operation
- Anomaly 1: Generator fault
- Anomaly 2: Transmission line issue
- Anomaly 3: Load anomaly

Benefit: Prevent cascading failures

Parameters

Performance Metrics

Activation Keywords

• anomaly detection
• log-likelihood ratio
• dynamical systems
• industrial control
• fault detection
• multiple model estimation
• hypothesis testing
• error rate characterization

Related Skills

• automated-cps-testing-act: CPS testing framework
• data-driven-moving-horizon-estimation: State estimation
• koopman-representation-learning: Dynamical system analysis

References

• Paper: arXiv:2604.11631 (April 2026)
• Authors: Riveiros, Barreau, Bastianello
• Theory: Log-likelihood ratio, Kalman filtering
• Application: Industrial control systems

Example Usage

"Implement LLR-based anomaly detection for control systems"
"Characterize error rates for multi-model fault detection"
"Design observer for industrial process monitoring"
"Detect subtle anomalies in linear Gaussian systems"

Code Template

import numpy as np
from scipy.stats import multivariate_normal
from filterpy.kalman import KalmanFilter

class LLRAnomalyDetector:
    def __init__(self, models):
        """
        Initialize with dictionary of system models.
        
        Args:
            models: Dict of {'name': LinearSystemModel}
        """
        self.models = models
        self.filters = {
            name: self._build_kalman_filter(model)
            for name, model in models.items()
        }
        self.llr_accumulators = {name: 0 for name in models}
    
    def _build_kalman_filter(self, model):
        """Build Kalman filter from system model."""
        kf = KalmanFilter(dim_x=model.n_states, dim_z=model.n_outputs)
        kf.F = model.A
        kf.B = model.B
        kf.H = model.C
        kf.Q = model.Q
        kf.R = model.R
        return kf
    
    def update(self, observation, control):
        """
        Process new observation and detect anomalies.
        
        Returns:
            detected_model: Most likely model name
            llrs: Log-likelihood ratios for all models
            confidence: Detection confidence
        """
        llrs = {}
        
        for name, kf in self.filters.items():
            # Prediction step
            kf.predict(u=control)
            
            # Compute innovation
            innovation = observation - kf.H @ kf.x_prior
            innovation_cov = kf.H @ kf.P_prior @ kf.H.T + kf.R
            
            # Log-likelihood of observation
            ll = multivariate_normal.logpdf(
                innovation, 
                mean=np.zeros_like(innovation),
                cov=innovation_cov
            )
            
            # Update accumulator
            self.llr_accumulators[name] += ll
            llrs[name] = self.llr_accumulators[name]
            
            # Update Kalman filter
            kf.update(observation)
        
        # Detect most likely model
        detected_model = max(llrs, key=llrs.get)
        
        # Compute confidence
        llr_values = np.array(list(llrs.values()))
        confidence = np.max(llr_values) - np.mean(llr_values)
        
        return detected_model, llrs, confidence
    
    def reset(self):
        """Reset LLR accumulators."""
        self.llr_accumulators = {name: 0 for name in self.models}

Notes

• Requires accurate system models for each anomaly type
• Linear Gaussian assumption enables theoretical analysis
• Real-time capable with Kalman filtering
• Error rate characterization aids observer design
• Multiple models enable anomaly isolation