Weighted Information Filtering, Smoothing, and Out-of-Sequence Measurement Processing

We consider the problem of state estimation in dynamical systems and propose a different mechanism for handling unmodeled system uncertainties. Instead of injecting random process noise, we assign different weights to measurements so that more recent measurements are assigned more weight. A specific choice of exponentially decaying weight function results in an algorithm with essentially the same recursive structure as the Kalman filter. It differs, however, in the manner in which old and new data are combined. While in the classical KF, the uncertainty associated with the previous estimate is inflated by adding the process noise covariance, in the present case, the uncertainty inflation is done by multiplying the previous covariance matrix by an exponential factor. This difference allows us to solve a larger variety of problems using essentially the same algorithm. We thus propose a unified and optimal, in the least-squares sense, method for filtering, prediction, smoothing and general out-of-sequence updates. All of these tasks require different Kalman-like algorithms when addressed in the classical manner.


I. INTRODUCTION
In the standard setup of the Kalman filter (KF) an unknown state sequence is modeled as a stochastic process evolving through a linear system and driven by an external random process noise. The sequence is observed by a linear measurement channel contaminated by random measurement noise. The resulting estimator is an optimal, in the minimum mean-square error (MMSE) sense, recursive algorithm for the estimation of the state using the available measurements.
The process noise sequence may be interpreted in two ways. From a mathematical point of view, it is the actual generator of the state sequence. In this case, the state and measurement models are assumed to be known perfectly and state estimation becomes a mathematical problem of Bayesian inference of a hidden sequence from noisy measurements. From a practical perspective, on the other hand, process noise is the designer's tool to cope with system uncertainties. For example, a (possibly) nonlinear dynamics governing the actual system's behavior may be modeled using a simple, linear motion model with process noise having sufficiently high variance. Clearly, the closer the system model to the true dynamics, the better the filter performance is expected to be. However, since system uncertainties are inevitable, using process noise cannot be avoided. The process noise is, thus, a method for handling unmodeled uncertainties in the dynamical (and, in fact, measurement) models.
Numerous generalizations and modifications to the basic version of the KF have been proposed. Optimal prediction and smoothing algorithms [1] are only two examples. A more recent notion of out-of-sequence measurement (OOSM) processing has become popular since the publication [2]. In this setup, some measurements are delayed and arrive after other, more recent measurements have already been processed. It is then required to update state estimates using the older measurements without reprocessing the entire measurement sequence. In some cases, such as Kalman smoothing and onestep OOSM processing [2], the derived algorithms remain optimal. There are scenarios, however, in which preserving optimality without increasing computational requirements is no longer possible. The problem of multi-step OOSM processing is a good example of such a scenario. Typically, the difficulty in deriving an algorithm or maintaining its optimality are rooted in the correlations generated by the process noise sequence and earlier measurements. On the other hand, the process noise is an inevitable mechanism in order to treat mismatches in the dynamical model of the state sequence. Consequently, some existing solutions to the multistep OOSM processing problem, e.g., [3], [4], make approximations thus sacrificing optimality. Other methods that address this problem increase the computation burden by augmenting the state vector (see [5], [6]).
Despite their similar form, the resulting optimal (or suboptimal) algorithms are different in all of the aforementioned cases. For example, the gain matrix computed for the fixedlag smoothing problem is computed differently than the gain for utilizing one-step OOSM. In other words, it is required to use a unique algorithm for each of the above problems.
In this paper we propose a different mechanism for handling unmodeled system uncertainties. Instead of injecting random process noise, we assign different weights to different measurements so that more recent measurements are assigned more weight. In this setup, the state at any time is a deterministic function of the state at any other time. We formulate a weighted least-squares cost for the estimation of the state using the available measurements and derive the corresponding solution. In this regard, the proposed approach is a generalization of the Savitzky-Golay filter [7] for polynomial smoothing of noisy data. A specific choice of an exponentially decaying weight function results in an algorithm with essentially the same recursive structure as the Kalman filter. An important difference, however, is the manner in which old data are combined with newly obtained measurements. While in the classical KF, the uncertainty associated with the previous estimate is inflated in an additive manner (by adding the process noisy covariance to the state error covariance matrix), in our case, the uncertainty inflation is done by multiplying the previous covariance matrix by an exponential factor. This simple difference allows us to solve a larger variety of problems using essentially the same algorithm. These problems include filtering, prediction, smoothing and general out-of-sequence updates, all of which, as mentioned, require a unique Kalman-like algorithm when addressed in a classical manner. In addition, the proposed algorithm does not require additional computational effort when applied to more complex problems.
The remainder of the paper is organized as follows. In Section II we formally state the problem at hand and list the relevant assumptions. Section III is devoted to the derivation of several algorithms starting with the simplest case of a batch-type estimator and proceeding to more complex, recursive routines. In Section IV several important properties of the derived estimator are discussed and some illustrative simple examples are provided. In Section V we consider an important special case of an exponentially decaying weight function and compare the resulting algorithm with the classical KF. In Section VI several applications are considered, the main of which is OOSM processing in the filtering and smoothing framework. Concluding remarks are given in Section VII. Due to space limitations, the proofs of the theorems are omitted.

II. PROBLEM SETUP
We consider a dynamical system evolving in an arbitrary manner in time. The evolution of the system is modeled by the following deterministic equation where x k is the system state at time t k , A k,ℓ is the known state transition matrix from time t ℓ to time t k and u k,ℓ is the cumulative effect of the deterministic inputs applied during the time interval (t ℓ , t k ]. The initial state x 0 is assumed to be deterministic but unknown. Clearly, the simple modeling may deviate significantly from the true system dynamics. In this case, there will be a model mismatch. We emphasize that we do not assume that the system is driven by a (random) process noise meaning that the only uncertainty in the described system is the initial state x 0 . In other words, no probabilistic modeling is associated with the sequence {x k }. For all k, ℓ ∈ N, the matrices {A k,ℓ } are assumed to be invertible [2]. Similarly to (1), assuming appropriate controllability conditions, we may express x ℓ using x k as follows where A ℓ,k = A −1 k,ℓ and C ℓ,k u ℓ,k = −A −1 k,ℓ C k,ℓ u k,ℓ . This allows the formulation of a direct observation of x ℓ as follows Here H ℓ is a known matrix, and {v ℓ } is a sequence of independent Gaussian random variables such that v ℓ ∼ N (0, R ℓ ). Note that, using (2), (3) also serves as an indirect observation of x k . Our goal is estimating the state x k using the available measurements {y ℓ , ℓ ∈ L} for some set of indices L. Note that since the state sequence is deterministic, every measurement in the set L may be viewed as a measurement of x k without adding uncertainty (in addition to the measurement noise v ℓ ). This is different from the Bayesian setup in which, in addition to the measurement noise, there is additional uncertainty due to the contribution of the process noise in the time interval defined by the considered measurement and the desired state x k . In the absence of a probabilistic prior for x k , and taking into account the symmetry of the distribution of v k , we utilize the relation (2) and consider the following weighted least-squares (WLS) cost (4) In other words, (4) is a sum of normalized squared differences between the measurements y ℓ and linear predictions thereof given by The weighting matrix of the ℓ-th element in the above equation comprises two elements -the inverse covariance matrix of the corresponding measurement noise and an additional matrix W ℓ,k which will play a crucial role in the sequel. Intuitively, the dependence of W ℓ,k on k causes this matrix to account for the relevance of the distant measurement y ℓ to the estimation of the present state x k thus controlling the contribution of y ℓ in case of a model mismatch. Perhaps the simplest example of W ℓ,k is the following, scalar-valued Savitzky-Golay [7] kernel Additional discussion on W ℓ,k and its properties is presented in Sections III and V. Note that the cost (4) could be obtained for the following formulation of the measurement equation for some index k, whereW ℓ,k = W −T ℓ,k and the effective covariance is We note in passing that, although defining different from (3) measurements, this model yields a standard WLS cost with optimal weights meaning its solution is the best linear unbiased estimator.
In the following section we derive optimal, in the sense of (4), estimator and consider special cases in which recursive formulation of the estimator is possible.

A. Batch-Form Estimator
We begin with the derivation of the optimal estimator, in batch form, for the quadratic cost (4). For clarity of the exposition we consider (6) as the measurement model. The first result, which is a straightforward generalization of the classical LS solution, is summarized in the following two Theorems.
Theorem 1: Let {y ℓ , ℓ ∈ L} be a general set of measurements defined in (6). The optimal solution of the problem (4) subject to the constraint (1) iŝ (8) Theorem 2: The covariance of the estimator (8) is given by Equation (2) allows us to rewrite (8) in the following compact form For notational convenience and without loss of generality, in the sequel, we set u ℓ,k = 0 for all k ∈ N and ℓ ∈ L.

B. Recursive Estimator
We now consider a narrower family of the weighting matrices W ℓ,k that depend explicitly both on the time of the currently considered measurement y ℓ , and on the time of the state to be estimated x k , thus accounting for the relevance of y ℓ on the estimation of x k . As we show in the sequel, this family allows recursive computation of the estimator (8) and the corresponding covariance (9).
Recall the considered measurement model (6). The weight assigned to y ℓ reduces as the difference t k − t ℓ increases. In this work we take W ℓ,k to be a scalar-valued function with the following properties: The simplest, trivial example of a function satisfying the properties (11a)-(11c) is W ℓ,k ≡ 1. A non trivial example of such a function, to be used in the sequel, is Informally, the parameter τ may be viewed as a "weight decay factor" or as a "decorrelation interval" so that the amount of information about the state, stored in the measurements, reduces as the temporal distance between these measurements and the state increases. The existence and utility of additional forms of W ℓ,k will be considered in a future work.
The three main results of this paper are stated next. Theorem 3: Let L = {1, 2, . . . , k}. The estimator (8) with the weighting matrix W ℓ,k satisfying (11) may be computed recursively as followŝ where the estimator covariance, P k , is given by The setup of Theorem 3 implies that the measurements {y ℓ , ℓ ∈ L} arrive sequentially, such that y k is the measurement arriving at time t k thus being an observation of x k according to (3). We show next that this does not have to be the case. To this end, consider a more general set of indices such that the measurement y ℓ k arrives at time t k , but corresponds to time t ℓ k . Clearly, if ℓ k = k for all ℓ k ∈ L, the present setup reduces to the one discussed in Theorem 3. Theorem 4: Let L be as defined in (15) and let d denote the k-th index in L such that t d = t ℓ k . The estimator (8) with the weighting matrix W ℓ,k satisfying (11) may be computed recursively as followŝ where the estimator covariance, P k , is given by (17) It is important to emphasize the implications of Theorem 4. The optimal update of the estimator is not affected by the times at which the measurements are acquired, as long as the actual validity time of these measurements is known. This means that seemingly different problems such as filtering and out-of-sequence measurement processing may be treated in essentially the same manner using (16). In other words, measurements arriving out-of-sequence are optimally processed using the same algorithm as those arriving in-sequence. This is not the case in the Bayesian setup, where a different algorithm is required to incorporate OOSM. Moreover, except for some special cases, standard OOSM processing is either suboptimal or requires increased computational resources. On the contrary, in the discussed setup, OOSM processing remains optimal in the sense of Theorem 4.
Inspecting the filter equations of Theorem 4 we notice that additional generalization is in place. The first summand in both (16) and (17) stands for the contribution of the previously updated estimate to the current time. However, the time index k does not have to refer to a specific time in the past and may, in fact, refer to any time in the past, present, or future. Therefore, we may rewrite the equations of the presented algorithm so that, in addition to the "measurement update step" defined by the second summand in (16) and (17) that allows arbitrary measurement timing, the "time propagation step" is also arbitrary as summarized in the following Theorem.
Theorem 5: Let L be as defined in (15) and let d denote the k-th index in L such that t d = t ℓ k . In addition, let m refer to an arbitrary index representing time instance t m at whicĥ x m and P m are known. The estimator (8) with the weighting matrix (11) may be computed recursively as followŝ where the estimator covariance, P k , is given by Clearly, Theorem 5 generalizes both Theorem 3 and Theorem 4 which follow as special cases. At this point, we may further generalize the discussion that follows Theorem 4 by realizing that the problem of fixed-lag smoothing may be treated by the same algorithm as filtering and OOSM processing. While in the Bayesian setup this problem requires a modified algorithm (known as the fixed-lag smoother) and cannot be done using a standard Kalman filter, in our case the problems of filtering, fixed-lag smoothing and OOSM processing may be solved optimally using the same algorithm described in Theorem 5.

IV. PROPERTIES AND EXAMPLES
A. Examples 1) Static Parameter: Let A k,ℓ = I meaning that x k = x k−1 = . . . = x 0 . In this case (4) reads arg min We assign equal weights to all measurements in the set L by setting, for all ℓ, k, W ℓ,k = 1 and obtain the following form of the estimator (8): where the covariance, P k , is given by Note that W ℓ,k = 1 satisfies the set of properties (11). The corresponding recursive form of the estimator thus readŝ 2) Constant Acceleration Model: Consider a state comprising position, velocity and acceleration such that Here ∆t t k − t k−1 and H k = [1 0 0]. Assuming the measurement noise variance is constant and using W ℓ,k = 1 as before, (8) reduces tô which is the exact result of Savitzky-Golay [7] for polynomial smoothing of noisy data. This observation allows us to relate the Savitzki-Golay filter to the more general KF framework. In other words, polynomial smoothing can be formulated within a state-space model using appropriate transition matrix (1), and measurement geometry (3).

B. Properties
We now state several properties satisfied by the batch and recursive versions of the estimator.
Theorem 6: The estimator (8) is unbiased and the corresponding error covariance (9) attains the Cramér-Rao lower bound with equality rendering the estimator efficient.
Lemma 1: The optimal estimator in either batch (8), or recursive (13) form does not depend on the scaling of the measurement noise covariance. Lemma 1 states that knowing the measurement noise variance is only required up to a multiplicative constant. Thus, for a constant measurement noise variance, one does not need to tune the corresponding parameter in the filter. This desired property is only possible due to the lack of the process noise and does not happen in, e.g., Kalman-filter based algorithms.

A. Concept
In this section we consider a specific function W ℓ,k and discuss the resulting estimator. The function under consideration is (12) which is repeated here for convenience.
It is readily seen that this weighting function satisfies the set of properties (11) such that a recursive form of the estimator (8) is given by (13). From this point on, however, we drop the formulation (6) and return to the original setup (3). Note that this does not affect the optimality of the solution in the weighted least-squares sense (4). For this choice of W ℓ,k the recursive estimator (13) takes the following form where P k is given by (27) Hereafter, we term the above estimator Exponentially Weighted Information Filter (EWIF). Note that once we have dropped the formulation (6), P k given in (27) is no longer the covariance of the estimator, but, rather, an approximation. This is similar to a KF setup in which the model is not known precisely and an artificial process noise is used to compensate this mismatch, such that the computed covariance is not the true error covariance.
We note in passing that although exponential weighting has been addressed in the literature in the past [8], [1], this was never done as an approach to deterministically cope with model uncertainties. Thus, the setup considered in the present section is a novel approach to address state estimation problems under possible model mismatch.

B. Relation to Kalman Filter
In order to gain additional insight into the differences between the EWIF approach and a standard KF setup, we recall the standard KF equations.
where P k is given by the Kalman gain, K k reads and It will be convenient to consider the information form of the latter, which we express in a somewhat nonstandard form as followŝ where P k is given by Showing the equivalence of the two formulations is a matter of straightforward algebraical manipulations. Comparing (32) and (33) with, respectively, (26) and (27) we make the following observations. It is readily seen that if τ → ∞ in (26) and (27), e − t k −t k−1 τ → 1 and the dependence of x k on the measurement history becomes uniform so that no extra weight is assigned to the more recent measurements. In this regime the algorithm reduces to a standard KF by setting Q k−1,k = 0 in (32) and (33). On the contrary, if τ → 0, then e − t k −t k−1 τ → 0, andx k depends solely on the current measurement y k . This is similar to setting a very large process noise in the standard KF setup.
Further, in both cases, the algorithms linearly combine two types of data -historical information, captured by the previously obtained estimatex k−1 , and new information stored in the currently acquired measurement y k . In the KF case, the weight assigned to the history is controlled by the quantity (P k−1 + A k−1,k Q k−1,k A T k−1,k ) −1 which is the information (inverse of the covariance) associated with older measurements. The reduction in the information, relatively to P −1 k−1 from step k − 1, is dictated by the process noise covariance Q k−1,k such that in the case of a very unpredictable state sequence (large process noise), much uncertainty is added between two consecutive measurement updates and the amount of information stored in the previous measurements reduces. On the contrary, at the absence of process noise, the weight assigned to the history of the measurements in the EWIF approach is controlled by the quantity e − t k −t k−1 the information associated with the measurement history reduces similarly to the KF case. In other words, the difference between the KF-based approach and the proposed EWIF method boils down to how the contribution of historical measurements is managed -additive increase of the error covariance in the case of KF versus multiplicative inflation of the covariance in EWIF. Our final comment relates to specific information reduction mechanism. While in the standard KF setup the weight assigned to the measurement history is controlled by the state process noise or, equivalently, by our assumptions on the system predictability, in the proposed approach this weight control mechanism is completely decoupled from the state and only considers the measurement sequence. This exposes additional robustness of EWIF to system model mismatch, similarly to the Savitzky-Golay filter. The parameter τ becomes then nothing but a tuning parameter as opposed to the Bayesian case where the process noise covariance has to be carefully chosen for each state dynamics based on some prior knowledge.

C. Toy Example
We now demonstrate the utility of the EWIF approach in a simple synthetic example. We consider a one-dimensional sinusoidal signal observed by a sequence of linear measurements corrupted by zero-mean Gaussian noise with variance 0.5. Note that the state at k+1 cannot be described as a linear transformation of the state at k. Thus, in order to track the sequence with a linear filter one has to compensate the model mismatch. We estimate the sequence of states using a linear Kalman filter in which the model mismatch is compensated by a process noise such that the state equation is where {w k } is a zero-mean white process noise sequence with standard deviation σ w = 0.1. Alternatively, the states are tracked using the EWIF algorithm with decorrelation constant τ = 0.5. The resulting estimates are presented in Fig. 1. It is readily seen that, for the described parameter setup, both algorithms achieve similar performance. This demonstrates the history handling mechanism carried out by the two methods as discussed in Section V-B and also suggests that, in similar estimation problems, for a given value of σ w , one can find an "equivalent" value of τ .

VI. APPLICATIONS
In this section we consider several classical applications commonly addressed and show how these can be treated with the considered deterministic framework.

A. Prediction and Smoothing
In order to derive optimal, in the sense of Section III, predictor we consider (18) and (19) for some index that m representing t m at whichx m and P m are known and k is an index representing some future time t k t m + ∆. It is assumed that no additional informative measurements are available for the computation ofx k . This is equivalent to setting infinite covariance to any new measurement thus nullifying the second summand in (18) and (19) which then becomex Note that the state prediction equation is exactly the same as in the KF case and the difference in the covariance expression is along the lines of the discussion in Section V-B. Some smoothing problems are trivially obtained from (18) and (19). For example, fixed-lag smoothing is obtained similarly by considering a constant lag t k − t m in the above equations.

B. Out-of-Sequence Measurement Processing
In centralized multisensor fusion systems it is common practice that some measurements arrive "out-of-sequence", namely, after other measurements, with later time-tags, have already been processed. The main challenge in utilizing OOSM, in the Bayesian setup, is the correlation generated by the past measurements and the process noise sequence and significant effort has been done in this direction. Some recent contributions, in addition to those cited in Section I, include [5], [9] and references therein. In all cases, a significant modification of the KF is required in order to cope with OOSM either exactly or approximately. In [5], e.g., state augmentation is required, thus increasing the computation burden when utilizing OOSM.
Next, we demonstrate how the OOSM utilization may be done in a simple manner using the EWIF approach. First, we rewrite (16) and (17) using the weighting function (25).
where the estimator covariance, P k , is given by Recall that t k is the time at which the state update is required and the current measurement is acquired, t k−1 is the previous update time, and t d = t ℓ k is the time for which the currently acquired measurement, y ℓ k ≡ y d , is valid. If y ℓ k is an insequence measurement, i.e., t d = t k , then e − t k −t d τ = 1 and the update rule reduces to that of (13). If t d = t k , the contribution of y ℓ k is reduced by an exponential factor accounting for the temporal distance of y ℓ k from the current time. It follows that processing an OOSM differs from utilizing an in-sequence measurement by the computation of an exponent of a scalar and, practically, requires no extra computational effort. This does not occur in any Bayesian approach.
One may wonder whether the performance of the algorithm depends on the order the measurements arrive. Invariance to the order of data processing is a desirable property of an algorithm. It is clear that, in the Bayesian setup, MMSE-optimal processing (which is not always feasible) is invariant to the order data are processed. Observing (36) and (37) it becomes evident that the same occurs in the proposed approach. Indeed, both of the above equations depend solely on the measurement validity time and not on the time these measurements arrived. Another way to arrive at the same conclusion is by observing the batch form of the estimator (8). It is readily seen that different orders of data processing correspond to different orders the summations in (8) are performed. Clearly, all yield identical results.
In order to demonstrate the algorithm for OOSM, we consider an example in which target moves arbitrarily in 2D space. Its position is observed by two sensors, one of which is undergoing a delay. The noise variance and measurement frequency vary along the scenario. The delay durations change randomly. The same EWIF algorithm summarized in (36) and (37) is used to track the target under four measurement regimes -filtering and smoothing, with and without incorporating OOSM. A snapshot of the trajectory accompanied by the corresponding measurements and estimates is presented in Fig. 2. An animation of the complete scenario can be found at https://youtu.be/2G1nEoNwwjY.
A sample implementation of the EWIF algorithm for a simpler example of fixed-lag smoothing in the presence of OOSM may be found at https://tinyurl.com/y3e4df3s.   Fig. 2: Tracking a target with arbitrary dynamics under four measurement regimes using the same EWIF algorithm with τ = 2 (sec). A two-dimensional CA model is used as the target model. The blue thin line is the true trajectory. Pluses and circles are the in-sequence and out-of-sequence measurements, respectively. Dashed cyan and magenta curves are, respectively, the filter and 1-sec-lag smoother using only the in-sequence measurements. Black and green curves are, respectively, the filter and 1-sec-lag smoother using the insequence measurements and OOSM. An animation of the complete scenario can be found at https://youtu.be/2G1nEoNwwjY.

VII. CONCLUDING REMARKS
We presented a deterministic framework for state estimation in dynamical systems. As opposed to the classical approach, where (random) process noise is used to cope with inevitable uncertainties in the dynamical model, in the presented setup, this uncertainty is treated by assigning exponentially decaying weights to the observations, so that more recent measurements are assigned more weight than "older" ones. The state is then estimated using a WLS algorithm. For the chosen weighting function, the algorithm possesses a convenient recursive form similarly to the Kalman filter. Handling the inherent uncertainty in the system model by using the weighted information approach allowed us to treat several, seemingly different problems using essentially the same algorithm. These include filtering, fixed-lag smoothing and OOSM processing. Lastly, we refer to our first example to emphasize that it is not a-priori obvious whether the proposed approach is preferable over the classical usage of process noise and a definite answer is, probably, casedependent. A clear advantage of the present approach is its applicability to a variety of estimation problems resulting in the same algorithm.