Robust Singular Smoothers For Tracking Using Low-Fidelity Data

TL;DR

This paper introduces a unified robust smoothing framework that effectively handles low-fidelity, biased, and correlated data for tracking underwater autonomous platforms, improving analysis and potential real-time application.

Contribution

It combines singular covariance models with robust losses and constraints, addressing challenges in noisy, biased, and discretized data for navigation.

Findings

Robust modeling improves data analysis accuracy.

Framework performs well with low-frequency USBL data.

Potential for real-time tracking extension.

Abstract

Tracking underwater autonomous platforms is often difficult because of noisy, biased, and discretized input data. Classic filters and smoothers based on standard assumptions of Gaussian white noise break down when presented with any of these challenges. Robust models (such as the Huber loss) and constraints (e.g. maximum velocity) are used to attenuate these issues. Here, we consider robust smoothing with singular covariance, which covers bias and correlated noise, as well as many specific model types, such as those used in navigation. In particular, we show how to combine singular covariance models with robust losses and state-space constraints in a unified framework that can handle very low-fidelity data. A noisy, biased, and discretized navigation dataset from a submerged, low-cost inertial measurement unit (IMU) package, with ultra short baseline (USBL) data for ground truth, provides an opportunity to stress-test the proposed framework with promising results. We show how robust modeling elements improve our ability to analyze the data, and present batch processing results for 10 minutes of data with three different frequencies of available USBL position fixes (gaps of 30 seconds, 1 minute, and 2 minutes). The results suggest that the framework can be extended to real-time tracking using robust windowed estimation.

Tables

Table 1. TABLE I : Prox operators of common PLQ penalties.

Penalty $f$	${\mathrm{prox}}_{\alpha f}\left(z\right)$	Ref.
$\frac{1}{2}{\Vert x\Vert }^2$, Fig. 2(a)	$\frac{1}{1+\alpha }z$	[28, 29]
${\Vert x\Vert }_1$, Fig. 2(b)	$\text{sign}\left(z\right)\odot {\left(\left|z\right|-\alpha \right)}_{+}$	[30, 31]
$q_{\tau }$, Fig. 2(c)		[32, 33]
$h_{\kappa }$, Fig. 2(d)	$\frac{\alpha }{\alpha +\kappa }z+\frac{\kappa }{\alpha +\kappa }{\mathrm{prox}}_{(\alpha +\kappa )\parallel \cdot {\parallel }_1}\left(z\right)$	[34]
$q_{\tau ,\kappa }$, Fig. 2(e)	$\frac{\alpha }{\alpha +\kappa }z+\frac{\kappa }{\alpha +\kappa }{\mathrm{prox}}_{\left(\alpha +\kappa \right)q_{\tau }}\left(z\right)$	[35]
$v_{ϵ}$, Fig. 2(f)		[36]
hubnik-$\kappa$, Fig. 2(g)	$\frac{\alpha }{\alpha +\kappa }z+\frac{\kappa }{\alpha +\kappa }{\mathrm{prox}}_{\left(\alpha +\kappa \right)v_{ϵ}}\left(z\right)$	[37, 38]
e-net, Fig. 2(h)	${\mathrm{prox}}_{\frac{\alpha }{1+2\alpha }\parallel \cdot {\parallel }_1}\left(\frac{1}{1+2\alpha }z\right)$	[39]

Equations

$$x_1$$

$$x_k$$

$$y_k$$

$$x_k$$

$$\underset{x_1\in X_1,\dots ,x_n\in X_N}{\min }$$

$$+{\rho }_m(R_k^{-1/2}(H_k{x}_k-y_k)),$$

$$\underset{x,u,t}{\min }$$

$$Q_k^{1/2}{u}_k=G_k{x}_k-x_{k-1}$$

$$R_k^{1/2}{t}_k=H_k{x}_k-y_k$$

$${\rho }_s(x_k)=\{\begin{array}{l}0x_k\in X_k\\ \infty x_k\notin X_k.\end{array}$$

$$\underset{z}{\min }$$

$$\rho (z)$$

$$z^T$$

$${\hat{w}}^T$$

$$A=\left(\begin{array}{cccc}{D}_1& 0& \dots & 0\\ B_1& D_2& 0& ⋮\\ 0& \ddots & \ddots & 0\\ 0& 0& B_{N-1}& D_N\end{array}\right),$$

$$D_i=\left(\begin{array}{ccc}{Q}_i^{1/2}& 0& I\\ 0& R_i^{1/2}& H_i\end{array}\right),B_j=\left(\begin{array}{ccc}0& 0& -G_{j+1}\\ 0& 0& 0\end{array}\right).$$

$$f^{\ast }(y)=\underset{x}{\sup }⟨x,y⟩-f(x),$$

$${\mathrm{prox}}_{\alpha f}(\zeta )=\arg \underset{x}{\min }\frac{1}{2\alpha }\Vert \zeta -x{\Vert }^2+f(x).$$

$${\mathrm{prox}}_{\rho }(z)+{\mathrm{prox}}_{{\rho }^{\ast }}(z)=z.$$

$$g(z)=\{\begin{array}{l}0Az=\hat{w}\\ \infty Az\ne \hat{w}\end{array}$$

$$\underset{z}{\min }\rho (z)+g(z)$$

$${\mathrm{prox}}_g(\eta )=\underset{Az=\hat{w}}{\mathrm{argmin}}\frac{1}{2}||\eta -z|{|}^2$$

$$\left[\begin{array}{cc}I& A^T\\ 0& A{A}^T\end{array}\right]\left[\begin{array}{c}z\\ \nu \end{array}\right]=\left[\begin{array}{c}\eta \\ A\eta -\hat{w}\end{array}\right]$$

$$\left[\begin{array}{ccccccc}0& \dots & 0& \underset{i}{\underbrace{1}}& 0& \dots & 0\end{array}\right]$$

$${\tilde{x}}_k=\left[\begin{array}{c}{x}_k\\ b\end{array}\right],{\tilde{Q}}_k=\left[\begin{array}{cc}Q& 0\\ 0& 0\end{array}\right].$$

$$w_k=M{w}_{k-1}+{\beta }_k,{\beta }_k\sim N(0,Q).$$

$${\tilde{x}}_k=\left[\begin{array}{c}{x}_k\\ w_k\end{array}\right],{\tilde{G}}_k=\left[\begin{array}{cc}{G}_k& I\\ 0& M\end{array}\right],{\tilde{Q}}_k=\left[\begin{array}{cc}0& 0\\ 0& Q_k\end{array}\right].$$

$$x_s=[x,y,z,\dot{x},\dot{y},\dot{z},\ddot{x},\ddot{y},\ddot{z}{]}^{\top }.$$

$${\dot{x}}_s$$

$$x_{s_{k+1}}$$

$$F_{s_k}$$

$$\Gamma =\left[\begin{array}{ccc}\frac{1}{3!}I{T}^3& \frac{1}{2!}I{T}^2& IT\end{array}\right]$$

$$R(\phi )=R_h^{\top }{R}_p^{\top }{R}_r^{\top },$$

$$\left[\begin{array}[]{ccc}c{h}&s{h}&0\\ -s{h}&c{h}&0\\ 0&0&1\end{array}\right],\quad\left[\begin{array}[]{ccc}c{p}&0&-s{p}\\ 0&1&0\\ sp&0&cp\end{array}\right],\quad\left[\begin{array}[]{ccc}1&0&0\\ 0&c{r}&s{r}\\ 0&-s{r}&c{r}\end{array}\right]$$

$$H_{s_k}=\left[\begin{array}{cc}{I}_{3\times 3}& 0_{3\times 6}\\ 0_{3\times 6}& R(\phi )\end{array}\right],z_{s_k}={\left[\begin{array}{cccc}{\mathrm{usbl}}^{\top }& {\ddot{x}}_{\text{meas}}& {\ddot{y}}_{\text{meas}}& {\ddot{z}}_{\text{meas}}\end{array}\right]}^{\top }$$

Full text

I Introduction

State-space models are ubiquitous in signal processing, and allow integration of disparate measurements to inform estimation, decisions, and control. Classic filtering [1, 2] and smoothing [3, 4] are core tools used to estimate these models. Their dependence on high-fidelity data, driven by Gaussian assumptions on errors and innovations, has demanded unequivocal attention from researchers and practitioners, and inspired robust dynamic inference methods.

While early robust approaches [5, 6] sought to modify iterations of the Kalman filter (KF) and Rauch-Tung-Striebel (RTS) smoother, over the last 25 years researchers have used robust formulations to weave assumptions on errors and innovations directly into the estimation problems themselves [7, 8, 9, 10, 11, 12, 13]. Constraints, when available, are also readily incorporated into the problem formulation [14]. Specifying the formulation leaves one free to choose from a range of optimization algorithms; the survey [15] describes a general class of models as well as first- and second-order methods to solve them.

Our focus is on models with singular variances for process and measurement residuals. These models are excluded by the assumptions of generalized smoothing [15] and all of the various special cases cited in the survey. In this paper, we build on the recently proposed framework of [16] for singular models, and systematically develop complementary modeling elements: robust penalties, informative constraints, and singular models. The resulting approach exploits the structure of singular covariances head on rather than using workarounds such as pseudo-inverses or variance boosting that either do not work in the general setting or introduce unnecessary changes to the fundamental model (see discussion in [16]).

A simple synthetic tracking example with a singular process model shows how the common tack of replacing the inverse by the pseudo-inverse fails dramatically in the presence of outliers (Figure 1). Any robust approach requires control of the null spaces associated to process and observations. We develop a direct practical formulation and method, and test it on a batch smoothing analysis of real field data, with a view towards real-time implementation in future work.

Background. Our main goal is to infer an unobserved state sequence $x_{1},\dots,x_{N}$ from noisy observations $y_{1},\dots,y_{N}$ using the model:

[TABLE]

where $x_{0}$ is a given initial state estimate, $G_{k}$ and $H_{k}$ are linear process and measurement models, $y_{1},\dots,y_{N}$ are observations, and $X_{k}$ specify additional information through constraints. The framework of [15] assumes that $w_{k}$ and $v_{k}$ are mutually independent random variables with known nonsingular covariances $Q_{k}$ and $R_{k}$, and that they follow from log-concave distributions; in particular they may be non-Gaussian.

Synthesizing all of this information gives the problem

[TABLE]

with $\rho_{p}$ and $\rho_{m}$ convex penalties. Using (2) provides estimates that are robust to outliers and can follow sudden changes in the state. Most of the inference- or optimization-based work in the convex dynamic setting is a special case of (2). Many examples, including robust penalties and constraints, are collected in [15].

We now extend to singular covariances $R_{k}$ (for errors $v_{k}$) and $Q_{k}$ (for innovations $w_{k}$). These models specify key use cases, particularly for innovations (process) modeling, briefly summaried below (see [16] for a more detailed discussion).

Deterministic integrals. Most models in robotics, particularly in navigation, use integration to model process relationships between state variables (e.g. when position, velocity, and acceleration are part of the state). Any deterministic integral yields a singular process model. The simplest example (with position a direct integral of stochastic velocity) is used to create Figure 1.

Nuisance parameters. Unknown constants that need calibration (such as fixed instrument biases) require special modeling in the nonsingular paradigm (2). With singular models, we can augment the state in order to infer these parameters.

Auto-regressive models and correlated errors. State-space models are broadly used in auto-regressive, moving average, and time series models [17]. These elements also appear in general smoothing models, particularly to deal with correlated measurement errors [18].

All three examples are accessible in the classic linear Gaussian setting. The KF need not invert $Q$ or $R$, and provides the minimum variance estimate for both singular and nonsingular models [19]. Some algorithms rely on precise knowledge of the error structure or explicit equality constraints [20, 21, 22]. Correlated errors are dealt with by augmenting the state and using a singular model [18]. None of these techniques generalize to singular models in the setting of (2), and some naive generalizations fail dramatically (Figure 1).

This paper builds on the reformulation of [16] for singular models. We develop a systematic approach and test it using a navigation model with a real-world mooring dataset. We show how robust statistics, singular models, and constraints can be systematically used to overcome a range of challenges simultaneously present in the dataset: (1) outliers, (2) deterministic relationships between states, (3) measurement biases, and (4) coarsely discretized observations.

The paper proceeds as follows. Section II summarizes the singular formulation of [16] and relevant optimization algorithms. Section III develops the key modeling elements to address common data challenges. Section IV presents the navigation models. Section V shows how the model elements come together to analyze the target mooring dataset, and obtain a high fidelity track from low-fidelity observations. Section VI concludes with discussion and future work.

II Robust Singular Formulation and Algorithm

We reformulate the robust smoothing problem to seamlessly allow both singular and nonsingular covarinace models for errors and innovations. The resulting problem can be solved with any primal-dual algorithm. We show how to apply the classic Douglas-Rachford splitting (DRS) algorithm (see e.g. [23, 24]) to the reformulated problem.

Problem (2) can be reformulated by introducing auxiliary variables $u_{k}$ and $t_{k}$ to represent pre-whitened innovations and measurement residuals:

[TABLE]

where $\rho_{s}(x_{k})$ may be taken as the convex indicator function to recover the constraints in (2):

[TABLE]

When $Q_{k}$ and $R_{k}$ are invertible, we can solve for $u_{k},t_{k}$ and recover (2). Otherwise, problem (3) is well-posed while (2) is not. We can write (3) in compact form

[TABLE]

where

[TABLE]

[TABLE]

and

[TABLE]

The variables are ordered in such a way that $A$ is block bi-diagonal. If all observations $z_{i}$ lie in the range of $H_{i}$, the constraint $Az=\hat{w}$ will be feasible [16].

The problem (4) is a convex optimization problem and can be solved using a variety of techniques. We show that the DRS algorithm is straightforward to implement, and preserves the computational complexity of the classic KF/RTS algorithms because of the way $A$ is structured in (6).

Given a convex function $f$, its convex conjugate $f^{*}$ is given by

[TABLE]

and its proximal operator with step $\alpha$, denoted by $\operatorname{prox}_{\alpha f}$ (see e.g. [25]) is given by:

[TABLE]

For a long list of objectives, prox operators are available in closed form or are efficiently computable. In particular this is the case when $\rho_{p},\rho_{m},\rho_{s}$ form any subset of the numerous elements briefly surveyed in Section III. It is actually the prox of $\rho^{*}$ that appears in the DRS iteration (Algorithm 1) rather than the prox of $\rho$, but these are linked by the simple formula

[TABLE]

To specify the algorithm, we let $g(z)$ be the indicator of the affine feasible region $Az=\hat{w}$:

[TABLE]

Problem (4) can now be written simply as

[TABLE]

which is a natural template for DRS, detailed in Algorithm 1.

To implement Algorithm 1 we need proximal operators of $\rho_{p},\rho_{m}$, and $\rho_{s}$. Eight common piecewise linear-quadratic (PLQ) penalties are shown in Figure 2, and their proximal operators are summarized in Table I.

The proximal operator for $g$ is given by

[TABLE]

which is a least squares problem with affine constraints. Solving it efficiently leverages the structure of (6). In particular we need to solve a single structured linear system

[TABLE]

where $AA^{T}$ is block tridiagonal and does not change between iterations. In our implementation, we need only compute a single block bidiagonal factorization once, which can then be used to solve (8) in $O(n^{2}N)$ operations in each iteration, no more expensive than a single matrix-vector multiply.

For piecewise-linear quadratic $\rho$ [26, 15], DRS converges to an optimal solution at a local linear rate [16], which does not depend on the condition number of $A$. A good initialization makes DRS competitive with the fastest available solvers, even second order methods with quadratic local rates [15].

III Modeling Elements

The proposed framework has three complementary modeling elements: singular covariance matrices $Q$ and $R$; process/measurement penalties $\rho_{p},\rho_{m}$; and constraints $\rho_{s}$ on the state. In this section, we show a range of choices for each element, and compute the operators required for Algorithm 1.

Singular covariances can be used to capture affine constraints, auto-regressive structure, integrated errors, and bias.

• •

Affine constraints using singular $R$. the $i$th element of the state at time $k$ is known exactly, add row

[TABLE]

to the measurement model $H_{k}$, a row and column of zeros to $R_{k}$, and the known value as the last element of $z_{k}$.

• •

Bias with singular $Q$. A common model for bias is to include it as a non-varying component across the state:

[TABLE]

• •

Correlated noise using singular $Q$. Correlated noise $w_{k}$ is typically modeled by [27]

[TABLE]

Here too we can augment the state and use a singular process variance:

[TABLE]

Piecewise linear-quadratic (PLQ) Penalties. The proposed framework allows process innovations, measurement residuals, and state regularization to come from any convex prox-friendly penalty. To keep the exposition simple, we collect eight commonly used convex piecewise linear-quadratic penalties in Figure 2, and compute their prox operators in Table I. The penalties can be thought of in terms of three features:

• •

Behavior at origin: nonsmooth features encourage exact fitting of the quantity being measured, while deadzones are appropriate for discretized observations.

• •

Tail growth: asymptotically linear penalties are more tolerant of large inputs. Applied to measurements, this gives robustness to outliers; applied to innovations, it gives an ability to quickly track evolving trends.

• •

Asymmetry: allows handling of special cases where under-estimating is qualitatively different from over-estimating.

Constraints. It is very convenient to enforce simple constraints on the state estimates $x_{k}$. If we take $\rho_{s}(x)=\delta_{X}(x)$ then the prox operator $\operatorname{prox}_{\rho_{s}}$ is simply the projection onto the set $X$. Box constraints are a very common type of constraints that enforce known bounds on the state, and have a trivial projection. The proposed framework allows us to use any convex region that has a computationally efficient projection.

IV Navigation Model

We use a constant-velocity kinematic model that is appropriate for many underwater vehicle applications, where accelerations are heavily damped and trajectories are often long straight lines (e.g. for transit or survey work). When the attitude is known or changing slowly, the model can be linearized effectively. For a vehicle that is well-instrumented in attitude, the uncertainty in position (and the x-y states in particular) is typically orders of magnitude larger than the uncertainty in attitude. Thus, in practice, we often simplify the full nonlinear vehicle process model to track only position states $(x,y,z)$, while assuming that the attitude states $(r,p,h)$ are directly available from the most recent sensor measurements. To make the model linear, the position and its derivatives are referenced to the local-level frame.

An effective model must counteract biases, outliers, and data discretization in the IMU data. We develop this model using the elements of the proposed framework.

Process model. To incorporate linear acceleration measurements from an IMU, we track linear velocities and linear acceleration in the state vector:

[TABLE]

The linear kinematic process model is given by

[TABLE]

where $w_{s}\sim\mathcal{N}(0,Q_{s})$ is zero-mean Gaussian noise. The linear process model (16) is discretized using a Taylor series:

[TABLE]

where $I$ in (17) denotes the $3\times 3$ identity matrix, and the higher order terms are identically zero because of the structure of $F_{s}$. We model the process covariance as if the error were the next term in the Taylor series approximation, a technique suggested by [40]. More precisely we set covariance to be the outer product, $\Gamma^{T}\Gamma$ where

[TABLE]

This leads to a rank 3 covariance for a $9\times 9$ matrix for a model that comprises position, velocity, and acceleration.

Given this covariance structure, the process model will penalize changes in acceleration. As the vehicle travels in a relatively straight line with small corrections, we expect to see acceleration mostly constant with a few small jumps. We use the $\ell_{1}$ norm for innovations, as it encourages exact fits while simultaneously allowing occasional sudden changes.

Measurement model. The inertial measurement unit (IMU) measures linear and angular accelerations relative to the physical frame of the vehicle on which it is mounted, while the state tracks linear acceleration relative to the navigation frame. We obtain the coordinate transformation between these frames using heading, pitch, and roll of the vehicle:

[TABLE]

where $R_{h}$, $R_{p}$, and $R_{r}$ are given by

[TABLE]

with $c\cdot$ and $s\cdot$ shorthand for $\cos(\cdot)$ and $\sin(\cdot)$.

Position data from the USBL is sampled at a lower update rate than the IMU. For any $s_{k}$ where such position data is available, we have the measurement model

[TABLE]

If there is no position data measured at time $s$ then we use the model

[TABLE]

The covariance used for measurement data similarly depends on whether there is position data available:

[TABLE]

where $U$ is a diagonal matrix reflecting position uncertainty, while $r_{s}$ captures uncertainty in IMU measurements.

Bias model. To compensate for the bias in acceleration data we augment the state vector to include bias variables:

[TABLE]

where $b_{1},b_{2},b_{3}$ are bias terms for acceleration in the $x,y,z$ directions in the local frame. To pass the bias estimates forward in time the process matrix is augmented with an identity block.

[TABLE]

At the first time point ($s_{k}=1$) we augment the covariance matrix with an identity block, and at all other time points we augment with a zero block. This adds equality constraints for the bias terms over all time points. The approach generalized easily to model piecewise-constant biases over longer periods.

[TABLE]

The measurement matrices are augmented with an identity block that shifts the acceleration measurements using the bias.

Discretization model. The measurement loss function is chosen to account for the level of discretization present in the data. The Vapnik loss (Figure 2f) has a ‘deadzone’ around the origin where small discrepancies are not penalized. This region is set using the level of discretization in the data, in this case $0.05$. The sharp corners of the Vapnik loss encourage errors to lie exactly in them, an unnecessary artifact. Thus we use the ‘Huberized’ version dubbed ‘hubnik’ (Figure 2 g). The prox operators for both losses are computed in Table I.

V Analysis of Mooring Data

We are interested in the ability to maintain an accurate position estimate on-board an autonomous underwater vehicle in real-time using acceleration measurements from a low-cost IMU, given periodic position fixes.

To test this, we use the singular robust Kalman framework to analyze data collected from a surface mooring equipped with an IMU on the subsea node. The mooring node, which is drifting with the current, is used as a proxy for a slowly moving underwater vehicle subject to unknown disturbances. We look at the position uncertainty and error accrued over time between the periodic, world-referenced position fixes provided by the USBL system.

Data description. Position fixes are available from a ship-based Sonardyne Ranger 2 USBL every 2 seconds, which we subsample to varying degrees for the analysis. Linear acceleration data from an LSM303D 3-axis accelerometer was collected at $\sim 0.075$ m/s2 precision using a Raspberry Pi Zero.

A small snippet of vertical acceleration data, Figure 3, shows the relatively coarse discretization of the measurements along with a mean that is shifted away from zero indicating bias. Outliers also appear likely. Because the IMU was generally upright with small perturbations in pitch and roll, the discretizations in the z-axis of the instrument frame are still clearly visible when the data are rotated into the vertical world frame. The small perturbations in pitch and roll cause the slight variations visible within each discretization level.

To illustrate our model’s usefulness on a data set such as in Figure 3, we isolate a 25-second window of the data and add modeling elements one at a time, noting the improvements they provide. We finish by running the full model with varying amounts of infrequent position data on a 10-minute section of data as a more practical experiment.

Model elements applied to data. We consider 25 seconds of IMU data and apply a classic Kalman filter (using least squares measurement and process loss) with the singular navigation model detailed in section IV. In practice, one would expect an underwater vehicle to be well instrumented in depth, but for illustration purposes we focus on depth and the vertical acceleration measurements. These measurements are the most biased and therefore improvements made by adding modeling elements are most clearly shown.

The 25 seconds of IMU data starts at a USBL measurement. Position is initialized to this starting USBL measurement, but for ease of comparison that position is treated as (0,0,0) and all subsequent USBL measurements are treated as relative offsets. DRS works better when its reasonably initialized. The initial state vector over the smoothing window is populated by setting acceleration to [math], and then propagating forward the most recent position fix with a significantly damped measurement of the most recent velocity to prevent divergence of the initial vector.

We start by adding bias estimation (for the acceleration measurements). This requires a second USBL position fix, 19 seconds into the 25-second data series. Figure 4 shows position and velocity estimates from a classic Kalman smoother applied to depth acceleration data with and without additional USBL position data. As expected due to the quality of the data we see poor performance in both cases although the small amount of additional position data offers slight improvement.

Figure 5 shows acceleration estimates of the classic Kalman smoother with and without additional USBL data. The two USBL measurements affect acceleration estimates locally but are generally overpowered by the vastly more complete set of (biased) acceleration data.

The picture improves significantly when we add in bias estimation (via singular process models), using the same two observations (compare Figures 4 and  6). With the bias removed, we can now focus removing the noticeable outliers form the acceleration data. For this purpose, we use the robust hubnik (Figure 2g) as our measurement loss. The ‘deadzone’ is designed to work with the coarse discretization of the acceleration data.

Figure 6 compares the results of (debiased) fitting between the least squares Kalman smoother and the robust version using the hubnik loss. Figure 7 shows the acceleration estimates with bias removed in both the classic and robust setting. There are noticeable differences in the acceleration estimates at around 11 and 22 seconds. When using the robust smoother, the effect of the outliers on the model’s estimate is greatly reduced (see the velocity estimates at 11 and 22 seconds). In the context of real-time tracking or forecasting, these sudden jumps will yield inaccurate predictions.

10 minutes of data. We apply the robust smoother with with bias, outliers, and discretization modeling elements to 10 minutes of IMU data. USBL data is available approximately every $2$ seconds, but we test performance of the smoothing algorithm at larger gaps, with USBL data supplied at $30,60,$ and $120$ seconds.

Figure 8 has the fitted position plots for all three frequencies. We see that without any USBL fix data, the estimate suffers; but even infrequent fixes give significant improvements when using the full capability of the smoother. In fact there are diminishing returns in increasing the USBL frequency; this is a promising result towards the future goal of a practical online implementation, particularly in settings where high-quality USBL observations are unavailable (e.g. during dives).

Figure 9 shows the fitted velocity. Here the effect of additional USBL data are more apparent, as velocity is completely inferred from position and acceleration. However, the smoothing estimates of velocity at infrequent USBL fixes are still very good compared to those informed by frequent USBL fixes.

VI Discussion and Future Work

We have proposed a singular Kalman smoothing framework that can use singular covariance models for process and measurements, convex robust losses, and state-space constraints. The modeler can use any convex loss that has an implementable prox, a class that includes the most common choices used for inference in tracking and navigation. The framework offers a wide range of flexibility that can be used to either counteract undesirable characteristics present in data or be used to increase model performance based off of relevant field knowledge. Future work will consider real-time implementation, as well as extension to nonlinear models.

Numerical experiments show that these modeling elements yield significant improvement on a noisy, challenging dataset. We also see that having a robust model makes the smoother less reliant on frequent high-quality position updates, which is a very promising development for underwater navigation.

This paper develops several tools required to move to robust singular tracking in real-time. A promising aspect of singular noise models is that they make it possible to do simple robust windowed smoothing, where estimates are constrained between windows as the tracking proceeds. Constraints on the state may play a bigger role in real-time estimation, since they can help detect outliers faster. Finally, robust penalties that provide better estimates may further improve performance of the DRS algorithm, by providing an effective initialization for each new window. We will focus on these developments in future work.

Acknowledgements

This material is based upon work supported by the Defense Advanced Research Agency (DARPA) and Naval Information Warfare Center (NIWC) Pacific under Contract No. N66001-16-C-4001. It is approved for public release, distribution unlimited. The views, opinions and/or findings expressed are those of the authors and should not be interpreted as representing the official views or policies of the Department of Defense or the U.S. Government.