Understanding the Kalman Filter

The Kalman filter is a recursive algorithm that estimates the state of a dynamic system from a series of noisy measurements. It is widely used in signal processing, control systems, navigation, and econometrics because it provides optimal estimates in the presence of Gaussian noise. The filter operates in two steps: prediction—where the next state is projected forward using a system model—and update—where the predicted state is corrected with the latest measurement. This iterative process yields estimates that are both accurate and computationally efficient.

At its core, the Kalman filter assumes the system follows a linear stochastic difference equation. The state evolves according to a known transition model, and measurements are linear functions of the state corrupted by additive noise. The filter maintains the first two moments of the state distribution—the mean and the covariance—which are sufficient for Gaussian distributions. By propagating these moments through the prediction and update steps, the Kalman filter minimizes the mean squared error of the state estimate.

Mathematically, the standard Kalman filter equations are:

  • Prediction: \(\hat{x}_{k|k-1} = F \hat{x}_{k-1|k-1} + B u_k\), \(P_{k|k-1} = F P_{k-1|k-1} F^T + Q\)
  • Update: \(K_k = P_{k|k-1} H^T (H P_{k|k-1} H^T + R)^{-1}\), \(\hat{x}_{k|k} = \hat{x}_{k|k-1} + K_k (z_k - H \hat{x}_{k|k-1})\), \(P_{k|k} = (I - K_k H) P_{k|k-1}\)

Where \(F\) is the state transition matrix, \(H\) is the observation matrix, \(Q\) is process noise covariance, \(R\) is measurement noise covariance, and \(K\) is the Kalman gain. For a simple scalar system (as often used for signal denoising), these reduce to the one-dimensional form presented later.

System Modeling for Signal Denoising

To apply a Kalman filter to denoise a signal, you must define a state-space model that describes how the clean signal evolves and how noisy measurements are generated. For many practical signals (e.g., sensor readings, financial time series, or biological signals), a simple constant-velocity or random-walk model suffices. The choice of model directly affects filter performance.

The most common model for denoising is the random walk:

  • State equation: \(x_k = x_{k-1} + w_k\), with \(w_k \sim \mathcal{N}(0, Q)\)
  • Measurement equation: \(z_k = x_k + v_k\), with \(v_k \sim \mathcal{N}(0, R)\)

Here, the state \(x_k\) is the unknown true signal at time \(k\), and the measurement \(z_k\) is the observed noisy sample. The process noise \(w_k\) accounts for gradual changes in the signal, while measurement noise \(v_k\) represents sensor or acquisition noise. If the signal is expected to change rapidly, a larger \(Q\) is appropriate; if the signal is nearly constant, a small \(Q\) works better.

For signals with known dynamics (e.g., sinusoidal oscillations, damped systems), you can use a higher-order model. For instance, a constant-velocity model tracks both position and velocity, which helps when the signal has a predictable rate of change. In MATLAB, you can represent these models using matrices and incorporate them into a general Kalman filter implementation.

Step-by-Step Implementation in MATLAB

Setting Up the Environment

Ensure you have MATLAB installed with the Signal Processing Toolbox (optional, but helpful for visualization). No special toolboxes are required for the basic Kalman filter; you can write the algorithm from scratch. Create a new script and define parameters at the top for easy tuning.

Defining the State-Space Model

For a scalar random walk model, the state transition is simply the identity (since \(F=1\)) and there is no control input (\(u_k=0\)). The observation matrix \(H\) is also 1. In MATLAB, these become:

F = 1;  % state transition matrix
H = 1;  % observation matrix
B = 0;  % control matrix (none)
u = 0;  % control input

For a two-state model (position and velocity), you would use:

F = [1 dt; 0 1];  % dt is the time step
H = [1 0];        % we only measure position
B = [0; 0];
u = [0; 0];

Adjust based on your signal characteristics. For this tutorial, we stay with the scalar model.

Initializing Parameters

You need initial estimates for the state (\(\hat{x}_0\)) and the error covariance (\(P_0\)). A common approach is to set \(\hat{x}_0\) equal to the first measurement and \(P_0\) to a large number (reflecting high uncertainty). For process and measurement noise covariances \(Q\) and \(R\), you often estimate them from the data or via trial and error. In MATLAB:

x_est = z(1);       % initial state estimate = first measurement
P = 1;              % initial estimation error covariance
Q = 0.01;           % process noise covariance (tune)
R = 0.1;            % measurement noise covariance (tune)

Tuning \(Q\) and \(R\) is critical. A smaller \(Q\) relative to \(R\) makes the filter trust the model more and smooth aggressively; a larger \(Q\) makes it follow measurements closely. See the advanced section for more guidance.

Prediction and Update Steps

The core loop iterates over each measurement. Inside the loop, we perform the prediction and update as defined by the Kalman equations. For the scalar case, this simplifies to:

N = length(z);
x_denoised = zeros(size(z));
x_denoised(1) = x_est;

for k = 2:N
    % Prediction
    x_pred = F * x_est + B * u;   % since F=1, B=0: x_pred = x_est
    P_pred = F * P * F' + Q;      % = P + Q
    
    % Update
    K = P_pred * H' / (H * P_pred * H' + R);  % = P_pred / (P_pred + R)
    x_est = x_pred + K * (z(k) - H * x_pred);
    P = (eye(size(P)) - K * H) * P_pred;
    
    % Store denoised value
    x_denoised(k) = x_est;
end

Note that for a scalar case, `eye(size(P))` is just 1. The denominator in Kalman gain uses `/` for scalar division. If you are using matrices, use `\` or `inv` appropriately, but for scalar it’s fine.

Example: Denoising a Noisy Sine Wave

Let’s test the filter on a synthetic signal. Create a clean sine wave, add Gaussian noise, then run the Kalman filter and compare results.

% Generate clean signal
t = 0:0.01:10;                 % time vector
x_true = sin(2*pi*0.5*t);     % clean 0.5 Hz sine wave

% Add noise
rng(42);                       % for reproducibility
noise_std = 0.3;
z = x_true + noise_std * randn(size(t));

% Kalman filter as defined above (with Q and R tuned)
Q = 0.001;                     % small process noise because sine is smooth
R = noise_std^2;               % measurement noise variance
x_est = z(1);
P = 1;
x_denoised = zeros(size(z));
x_denoised(1) = x_est;

for k = 2:length(z)
    % Prediction
    x_pred = x_est;
    P_pred = P + Q;
    
    % Update
    K = P_pred / (P_pred + R);
    x_est = x_pred + K * (z(k) - x_pred);
    P = (1 - K) * P_pred;
    
    x_denoised(k) = x_est;
end

% Plot results
figure;
plot(t, x_true, 'b-', 'LineWidth', 1.5); hold on;
plot(t, z, 'r.', 'MarkerSize', 4);
plot(t, x_denoised, 'g-', 'LineWidth', 1.5);
legend('True Signal', 'Noisy Measurements', 'Kalman Filter Output');
xlabel('Time (s)');
ylabel('Amplitude');
title('Kalman Filter Denoising of a Sine Wave');

This example demonstrates that the Kalman filter can nearly recover the original sine wave from noisy observations. The green output lags slightly due to the predictive nature of the filter, but this lag can be reduced by increasing \(Q\) (trusting measurements more) at the cost of more noise.

Advanced Considerations

Tuning Q and R

Selecting appropriate noise covariances is the most challenging aspect of using a Kalman filter. If you have prior knowledge of sensor noise, set \(R\) equal to the measurement variance. You can estimate \(R\) from static data by computing the variance of the signal when no dynamics are expected. For process noise \(Q\), consider the expected rate of change. A common heuristic is to set \(Q\) proportional to the expected squared change per time step.

Another approach is optimization: run the filter over a grid of \(Q\) and \(R\) values and choose the combination that minimizes the mean squared error between the filtered output and a known ground truth (if available). In production systems, adaptive Kalman filters adjust \(Q\) and \(R\) online based on innovation statistics.

Handling Non-Linear Systems with the Extended Kalman Filter

If your signal or measurements involve non-linear relationships, the standard Kalman filter is inadequate. You can use the Extended Kalman Filter (EKF), which linearizes the non-linear functions around the current estimate using Jacobians. MATLAB offers the `ekf` functions in the Control System Toolbox or you can implement it manually. For instance, if the measurement equation is \(z_k = \sin(x_k) + v_k\), you would compute the Jacobian \(H = \cos(\hat{x}_{k|k-1})\). EKF is widely used in robotics (e.g., for localization) and non-linear signal processing.

Performance Optimization

For real-time applications, ensure the filter loop is computationally efficient. Pre-allocate arrays and avoid dynamic resizing. Use vectorized operations where possible (e.g., for offline processing, you can implement a forward-backward smoother). MATLAB’s built-in `kalman` function from the Control System Toolbox can also speed up development. For large datasets, consider using the Unscented Kalman Filter (UKF) which avoids explicit Jacobian computation and often performs better for highly non-linear systems.

Common Pitfalls and How to Avoid Them

  • Poor initialization: Setting initial state too far from true value can cause slow convergence. Use the first measurement as a starting point.
  • Unrealistic covariances: If \(P\) becomes negative or non-positive definite due to numerical errors, the filter diverges. Use square-root forms or ensure proper factorization (e.g., Joseph form update). In MATLAB, use `1e-6` as a lower bound for variances.
  • Model mismatch: If the true signal dynamics differ significantly from the assumed model (e.g., using a random walk on a periodic signal with high frequency), the filter will not perform well. Consider increasing model order or using an adaptive approach.
  • Ignoring process noise when state is constant: If \(Q=0\) and the state is truly static, the filter may become overconfident and ignore measurements after many steps. Always include a small \(Q\) to keep the filter responsive.
  • Measurement outliers: Standard Kalman filters are sensitive to outliers. Use robust versions (e.g., filter with modified update or apply outlier detection before updating).

Applications of Kalman Filters Beyond Denoising

While this article focuses on signal denoising, Kalman filters are ubiquitous in engineering and data science. Common applications include:

  • Navigation and GPS: Fusing inertial measurement units (IMU) with GPS data to estimate position, velocity, and attitude.
  • Finance: Estimating the hidden state of asset prices (e.g., stochastic volatility) from noisy observed prices.
  • Robotics: Simultaneous localization and mapping (SLAM) often uses EKF.
  • Econometrics: Estimating unobserved components in time series (e.g., trend-cycle decomposition).
  • Biomedical signal processing: Denoising ECG or EEG signals, tracking eye movement.

For more comprehensive resources, refer to the MATLAB documentation on Kalman filtering, the Wikipedia article for theoretical background, or a classic textbook like "Optimal State Estimation" by Dan Simon.

Conclusion

Implementing a Kalman filter in MATLAB for signal denoising is a straightforward process once you understand the underlying assumptions and tuning parameters. By defining a suitable state-space model, initializing state and covariances appropriately, and iterating through prediction-update steps, you can effectively reduce noise and recover clean signals. The scalar random walk model is a good starting point, but you can extend it to more complex dynamics using matrix formulations. Experiment with different values of \(Q\) and \(R\) to balance smoothness versus responsiveness. With practice, the Kalman filter becomes an indispensable tool in your signal processing toolbox.