Modular Nonlinear Least Squares Filtering

• The paper introduces modular nonlinear least squares filtering by decomposing nonlinear estimation into independent modules like process modeling, observation handling, and local linearization.
• It presents recursive Gauss–Newton and Bayesian kernel methods that leverage online updates and MAP approximations to efficiently address nonlinear state estimation.
• The work details subsystem-wise updates using covariance intersection techniques to prevent information double counting while maintaining robustness under communication and computational constraints.

Modular nonlinear least squares filtering denotes a family of recursive state-estimation methods in which nonlinear filtering is formulated through least-squares structure, but implemented through separable modules such as a process model, an observation model, local linearization or moment computation, recursive information accumulation, and a correction rule. Across the works considered here, the term encompasses at least three closely related constructions: Jacobian-based recursive Gauss–Newton filtering, Bayesian kernel adaptive filtering interpreted as approximate sequential MAP estimation, and subsystem-wise nonlinear least-squares updates for coupled systems with Covariance Intersection to avoid double counting. A common theme is that a nonlinear estimation problem is decomposed into components that can be modified independently while preserving an online recursive architecture (Nadjiasngar et al., 2011, Park et al., 2013, Zamani et al., 24 Jul 2025).

1. Conceptual basis and modular decomposition

In the Recursive Gauss-Newton Filter, the filtering problem is posed as nonlinear state estimation for a state-space model

$D X(t) = F(X(t)), \qquad Y(t) = G(X(t)) + v(t),$

with Gaussian measurement noise, and the method is stated to apply to four cases: linear process / linear observation, linear process / nonlinear observation, nonlinear process / linear observation, and nonlinear process / nonlinear observation (Nadjiasngar et al., 2011). The derivation decomposes naturally into a process model module, an observation model module, local linearization through the process Jacobian

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$

and the observation sensitivity matrix

$$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$

together with recursive information accumulation and a nonlinear correction solve.

In the Bayesian reinterpretation of kernel least mean squares, modularity is expressed differently. The model is separated into latent dynamics $p(w_{k+1}\mid w_k)$, observation model $p(y_k\mid x_k,w_k)$, and an approximate inference rule that retains only the MAP estimate and discards covariance. The paper’s main message is that KLMS can be understood as an approximate Bayesian filtering algorithm, and that this viewpoint places it within a broader class of modular nonlinear least-squares filtering methods (Park et al., 2013). The same decomposition permits separate modification of state dynamics and likelihood.

In the subsystem-localization setting, modularity refers to decomposition over state blocks. The system is split into independent subsystems with separate state and covariance estimates, while relative measurements can depend on multiple subsystems simultaneously. The proposed approach updates each subsystem independently, even when a relative measurement couples them, and does so without tracking the cross-covariance $P_{12}$ (Zamani et al., 24 Jul 2025). In that formulation, modularity is architectural: state representation, covariance bookkeeping, and relative-measurement processing are all decomposed by subsystem.

These three uses of “modular” are not identical, but they are compatible. This suggests that modular nonlinear least squares filtering is less a single algorithm than a design pattern in which nonlinear least-squares structure is preserved while one selectively changes model components, approximation components, or subsystem boundaries.

2. Recursive Gauss–Newton formulation

The Gauss-Newton filter is derived from weighted least squares and Newton-style local linearisation. A nominal trajectory $\bar{X}(t)$ is introduced and the state is written as

$X(t) = \bar{X}(t) + \delta X(t),$

with first-order perturbation dynamics

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$

Similarly, the observation perturbation is expressed as

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$

Stacking $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$0 observations yields a linear weighted least-squares problem

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$1

whose solution is

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$2

The state update is then

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$3

This is the paper’s Gauss-Newton filter viewpoint: estimate the correction $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$4 by solving a locally linear weighted least squares problem (Nadjiasngar et al., 2011).

The recursive form eliminates the need to store the full batch. With forgetting factor $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$5, the information matrix and gradient-like vector are defined by

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$6

so that

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$7

The recursive updates are

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$8

$$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$9

Only the current measurement and current linearization are required; past information is compressed into $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$0 and $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$1.

The Levenberg–Marquardt adaptation modifies the normal equations by

$$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$2

with gain ratio

$$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$3

Large $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$4 yields a steepest-descent-like step, while small $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$5 recovers the Gauss-Newton step. The paper states that the recursive filter is robust in simulation, that the LM augmentation improves convergence and avoids singularity, and that recursion gives a large reduction in computational memory relative to the batch GNF (Nadjiasngar et al., 2011).

3. Bayesian kernel filtering and approximate sequential MAP inference

For kernel least mean squares, the filter is written in feature-space form

$$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$6

where $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$7, $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$8 is the feature map, and for an RKHS kernel $$M(\bar{X}_n)=\left.\frac{\partial G(X_n)}{\partial X_n}\right|_{\bar{X}_n},$$9,

$p(w_{k+1}\mid w_k)$0

The stochastic-gradient LMS/KLMS update is

$p(w_{k+1}\mid w_k)$1

and the resulting prediction can be written as

$p(w_{k+1}\mid w_k)$2

The paper’s central claim is that LMS is closely related to Kalman filtering, and therefore KLMS can be interpreted as an approximate Bayesian filtering method (Park et al., 2013).

The associated nonstationary latent-state model introduces time-varying weights: $p(w_{k+1}\mid w_k)$3

$p(w_{k+1}\mid w_k)$4

If the posterior is kept Gaussian,

$p(w_{k+1}\mid w_k)$5

then conjugacy yields

$p(w_{k+1}\mid w_k)$6

The paper notes that this exact Gaussian recursion is too costly in RKHS because the covariance is high-dimensional and grows with data.

To preserve KLMS-like efficiency, the approximation

$p(w_{k+1}\mid w_k)$7

is used. With this approximation and the Gaussian diffusion prior, the MAP recursion becomes

$p(w_{k+1}\mid w_k)$8

For normalized kernels, $p(w_{k+1}\mid w_k)$9, this reduces to

$p(y_k\mid x_k,w_k)$0

The classic KLMS update therefore emerges from a Gaussian diffusion prior and Gaussian observation model followed by a MAP approximation that discards covariance (Park et al., 2013).

This formulation provides an explicit explanation for the usual constant-step-size tracking behavior. The frequentist analysis in the same work studies a drifting true parameter

$p(y_k\mid x_k,w_k)$1

with steady-state relation

$p(y_k\mid x_k,w_k)$2

In a stationary environment $p(y_k\mid x_k,w_k)$3, smaller $p(y_k\mid x_k,w_k)$4 reduces steady-state error; in a nonstationary environment, finite $p(y_k\mid x_k,w_k)$5 is required for tracking. The paper interprets this as evidence that KLMS behaves like a filter for drifting latent parameters.

4. Modular extensions: forgetting, non-Gaussian observations, and quadratic measurement updates

Within the Bayesian KLMS framework, the first major extension changes the latent dynamics. Replacing the random walk by

$p(y_k\mid x_k,w_k)$6

yields a forgetting KLMS, denoted fKLMS, with update

$p(y_k\mid x_k,w_k)$7

It can also be expressed as

$p(y_k\mid x_k,w_k)$8

so the influence of old samples decays geometrically. The paper notes that this forgetting is introduced by changing the prior mean dynamics, not by enlarging covariance as in Kalman filtering (Park et al., 2013).

The second extension changes the observation model. For count data,

$p(y_k\mid x_k,w_k)$9

and for binary labels,

$P_{12}$0

In both cases, the MAP solution lies in the direction of $P_{12}$1, so the update can be written as

$P_{12}$2

and the optimization reduces to a scalar strictly concave problem in $P_{12}$3. The paper states that the complexity remains $P_{12}$4 with only a small scalar optimization overhead, and more generally that the MAP-only Bayesian approximations preserve $P_{12}$5 time per iteration, whereas extended KRLS or exact Bayesian covariance tracking is typically quadratic per update (Park et al., 2013).

A distinct line of modular extension appears in the quadratic EKF/UKF work. There, the prediction architecture is left unchanged, but the measurement update is upgraded from the linear estimator

$P_{12}$6

to the quadratic family

$P_{12}$7

The optimal quadratic estimator is

$P_{12}$8

with augmented gain

$P_{12}$9

The state and covariance updates preserve the Kalman form

$\bar{X}(t)$0

The paper emphasizes that EKF and UKF remain linear estimators because the corrected state is affine in the measurement residual, whereas QEKF and QUKF use the residual square and therefore produce a parabolic dependence on the data (Servadio et al., 6 Jun 2025).

In the scalar nonlinear measurement example

$\bar{X}(t)$1

the QEKF collapses to the EKF under Gaussian assumptions because $\bar{X}(t)$2 when skewness is zero, while the QUKF computes the required cross-moments from sigma points and improves accuracy by around $\bar{X}(t)$3 compared with the UKF. In the Clohessy–Wiltshire relative navigation example with nonlinear angle measurements and discrete non-Gaussian noise, both QEKF and QUKF are reported to be consistent and to outperform EKF/UKF (Servadio et al., 6 Jun 2025).

5. Subsystem-wise modular filtering and Covariance Intersection

The robot-and-landmark formulation provides an explicit nonlinear least-squares treatment of coupled subsystems. Two subsystems,

$\bar{X}(t)$4

have local priors $\bar{X}(t)$5 and zero initial cross-covariance,

$\bar{X}(t)$6

Relative measurements are represented by a nonlinear error map

$\bar{X}(t)$7

with weighted least-squares term

$\bar{X}(t)$8

For updating $\bar{X}(t)$9, the local objective is

$X(t) = \bar{X}(t) + \delta X(t),$0

After minimizing over $X(t) = \bar{X}(t) + \delta X(t),$1, the effect of the relative measurement becomes

$X(t) = \bar{X}(t) + \delta X(t),$2

where

$X(t) = \bar{X}(t) + \delta X(t),$3

The paper interprets $X(t) = \bar{X}(t) + \delta X(t),$4 as measurement covariance inflated by the uncertainty of the other subsystem pushed through its Jacobian (Zamani et al., 24 Jul 2025).

Because the relative-measurement-derived estimate can depend on an estimate that has itself already been influenced by the current subsystem, ordinary least-squares/Kalman fusion can double count prior information. The proposed remedy is Covariance Intersection: $X(t) = \bar{X}(t) + \delta X(t),$5

$X(t) = \bar{X}(t) + \delta X(t),$6

with

$X(t) = \bar{X}(t) + \delta X(t),$7

The same paper gives the least-squares interpretation

$X(t) = \bar{X}(t) + \delta X(t),$8

which allows CI to be embedded directly into the nonlinear least-squares framework rather than treated as a separate heuristic (Zamani et al., 24 Jul 2025).

For subsystem $X(t) = \bar{X}(t) + \delta X(t),$9, the modular update becomes

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$0

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$1

An analogous formula holds for subsystem $D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$2, with $D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$3. The architecture is modular in the specific sense that each subsystem keeps only its own $D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$4, relative measurements are processed by one-step elimination, and no global joint covariance is maintained.

6. Robot–landmark localization, performance trade-offs, and limitations

The robot-landmark specialization is posed on $D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$5. The robot pose is

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$6

with unicycle dynamics

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$7

while the landmark is stationary,

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$8

Robot prediction uses noisy velocity measurements and the EKF covariance recursion

$D(\delta X(t)) = A(\bar{X}(t))\,\delta X(t).$9

with occasional direct pose measurements

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$0

The relative bearing unit vector is modeled as

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$1

and the nonlinear coupling cost is

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$2

The residual is zero when the relative landmark vector aligns with the measured bearing, and only the component perpendicular to $$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$3 contributes (Zamani et al., 24 Jul 2025).

The modular landmark update and robot update are derived by absorbing the uncertainty of the other subsystem into an effective covariance and then applying CI-based fusion. The paper benchmarks this method against a monolithic joint-state estimator with state

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$4

in a randomized Monte Carlo study with a $$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$5 square region, horizon $$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$6, twist measurements at every step, robot pose measurements every 3 steps, bearing measurements every 6 steps, and 20,000 trials. The five reported methods are Joint, FSafe, FKalman, Safe, and Kalman. Using final landmark error

$$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$7

as metric, the reported mean and standard deviation are:

• Joint: mean $$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$8, std $$\delta Y_n = Y_n - \bar{Y}_n = M(\bar{X}_n)\,\delta X_n + v_n.$$9
• FSafe: mean $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$00, std $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$01
• FKalman: mean $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$02, std $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$03
• Safe: mean $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$04, std $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$05
• Kalman: mean $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$06, std $$A(\bar{X}(t))=\left.\frac{\partial F(X(t))}{\partial X(t)}\right|_{\bar{X}(t)},$$07

The paper states that Joint has the best overall inlier performance but more sensitivity/outliers, that FSafe is nearly as good as Joint on average with lower variance and fewer outliers, and that Safe and Kalman degrade gracefully rather than catastrophically under reduced communication and reduced computation (Zamani et al., 24 Jul 2025).

Across the cited works, several limitations recur. The Recursive Gauss-Newton Filter is fundamentally local, requires a nominal trajectory close enough to the true trajectory for linearization to remain valid, and still incurs iterative solve cost at each step (Nadjiasngar et al., 2011). The Bayesian KLMS extensions retain only the MAP estimate and discard covariance, so they do not quantify uncertainty over the weights and are not fully optimal Bayesian inference (Park et al., 2013). The quadratic EKF/UKF requires higher-order moments, and the paper notes that handling these moments can be delicate, especially with measurement noise and vector measurements; under Gaussian assumptions, the quadratic term may vanish, as in the QEKF scalar example (Servadio et al., 6 Jun 2025). The modular CI-based subsystem method sacrifices some optimality relative to a full joint filter, adds conservatism through CI, and depends on communication level (Zamani et al., 24 Jul 2025).

Taken together, these works support a technically specific understanding of modular nonlinear least squares filtering: nonlinear estimation is cast as weighted least squares or approximate MMSE correction, then organized so that individual modules—dynamics, observation law, Jacobians or moments, recursive information state, damping rule, approximate posterior representation, or subsystem fusion policy—can be altered without discarding the recursive filtering architecture.