Unscented Kalman Filter Empowered by Bayesian Model Evidence for System Identification in Structural Dynamics

System identification is often limited to parameter identification, while model uncertainties are disregarded or accounted for by a fictitious process noise. However, modelling assumptions may have a large impact on system identification. For this reason, we propose to use an unscented Kalman filter (UKF) empowered by online Bayesian model evidence computation for the sake of system identification and model selection. This approach employs more than one model to track the state of the system and associates with each model a plausibility measure, updated whenever new measurements are available. The filter outcomes obtained for different models are then compared and a quantitative confidence value is associated with each of them. Only the system identification outcomes related to the model with the highest plausibility are considered. While the coupling of extended Kalman filters (EKFs) and Bayesian model evidence was already addressed, we modify the approach to exploit the most striking features of the UKF, namely, the ease of implementation and higher-order accuracy in the description of the evolution of the state mean and variance. A challenging identification problem related to structural dynamics is discussed to show the effectiveness of the proposed methodology.


Introduction
Kalman filters (KFs) are well-known tools for system identification. They work by applying a predictor phase, in which a suitable model is needed to predict the evolution of a dynamic system, and a correction phase, in which corrections to the prediction are applied by recursively processing system measurements [1].
In civil and mechanical engineering, different model classes, consisting of different parametrizations of the structure to be identified, can be formulated. They are built upon different levels of complexity in the description of the system mechanics and uncertainty in the formulation of the modelling assumptions. Emphasis is usually placed on improving the quality of the parameter estimate, especially whenever nonlinear dynamic systems are handled. With this goal, KF extensions such as the extended Kalman filter (EKF) or the unscented Kalman filter (UKF) have been introduced. On the contrary, model uncertainties are often disregarded or accounted for by a fictitious process noise. In this work, we propose a way to tackle this aspect by calculating a quantitative estimate, referred to as model evidence, measuring how much the model employed by the KF is plausible with respect to other possible parametrizations. While a similar estimate was discussed in [2] for the EKF, here, we develop a model evidence formula suited for the UKF to exploit its ease of implementation and higher-order accuracy in the description of the evolution of the state mean and variance.
The remainder of the contribution is organized as follows: In Section 2, first, the governing equations of a mechanical elasto-dynamic system are discussed; second, the related algorithm showing the application of the UKF for parameter estimation is reported and finally, the equations allowing for recursive model evidence calculation are presented. In Section 3, a case study featuring a shear building excited by real ground acceleration is discussed, showing how parameter identification outcomes are affected by different structural parametrizations and how model evidence can be used for the sake of model selection a posteriori. Conclusions are finally discussed in Section 4.

Elasto-Dynamic Problem
We focus on situations where the system dynamics is described by the finite element (FE) discretised version of a general elasto-dynamic problem. At time t k+1 , it reads: where M, C and K are the mass, damping and stiffness matrices, respectively; q,q,q ∈ R n×1 are the nodal displacements, velocities and accelerations, respectively; and f k+1 ∈ R n×1 is the external force vector, assumed to be known. Equation (1) is integrated in time by using the α-method [3], ruled by the parameters α m , α f and β. At each time step, the displacement field q k+1 is obtained by solving The modified matrix K * and the right hand side vector f * k+1 are computed by where Moreover, the mechanical system is assumed to be only partially observed. Accordingly, a Boolean matrix H ∈ R n o ×3n establishes the connection between the n o observed quantitiesŷ k+1 ∈ R n o and the kinematic fields, as follows:

Unscented Kalman Filter for Parameter Estimation
In this study, the ultimate goal of filtering is to estimate the unknown parameters θ ∈ R n p ×1 ruling the mechanical response of the structure to be identified, where typically C = C(θ) and K = K(θ). In [1,4], Kalman filtering techniques were successfully applied, even in the presence of nonlinearities due to damage evolution in the observed system, by solving a dual estimation problem and thereby adopting as state variables the model displacements and the unknown parameters governing the response of the mechanical domain. However, treating FE solutions characterized by a large number n of degrees of freedom (DOF) may result in an excessive computational burden when dealing with dual estimation. A possible solution consists of obtaining a reduced order model (ROM) representation of the mechanical domain and adopting as state variables, instead of the nodal kinematics, the ROM DOF [5,6]. This strategy has been explored in [7]. Here, we consider only θ as a state variable to avoid the computational burden connected to the combined use of UKF and DOF tracking when large FE models are addressed, despite the enhanced performance usually guaranteed by state tracking [1,8]. The following state-space representation is used: where the θ is driven by a random walk ruled by w k ∈ R n p ×1 , modelled as a white process noise w k ∼ N (0, Q), and the FE predicted outputŷ k+1 is related to the actual response of the structure by adding a measurement noise v k+1 ∈ R n y ×1 , modelled as white v k+1 ∼ N (0, R). The matrices Q and R are symmetric and positive defined. The time variation of θ, introduced by the random walk formulation, is fictitious. KFs attempt to propagate the mean and the covariance of the state variable vector through the state-space and the measurement update equations. Instead of propagating the probability density functions associated with the state variables, it is indeed preferable to draw a set of samples, deterministically propagate them and finally compute the mean and covariance of the state vector; this is especially profitable when the state-space and/or the measurement update equations are nonlinear. The UKF is based on this idea. The propagated vector collects a set of so-called sigma points (SPs) ϑ i k , i = 1, . . . , (2n θ + 1), distributed such that the mean and covariance of these points match those of the state variables. A scaled version of the UKF is used by setting the parameters α SP , κ SP and β SP , as detailed in [9], to avoid sampling nonlocal effects that would spoil the state variable mean and covariance reconstruction [10]. In the predictor phase, this vector is propagated from the k-th to the k+1-th time step through the state-space equations. In the corrector phase, the estimated output covarianceP yy k+1|k and the estimated cross covarianceP θy k+1|k are used to compute the Kalman gain G k+1 needed to correct the propagated meanθ k+1|k and covarianceP θθ k+1|k on the basis of the collected measurements y k+1 . The full expression of these quantities and the application of the UKF are detailed in Algorithm 1, adapted from [11].

Model Evidence Computation for Unscented Kalman Filter
System identification is usually limited to select a particular parametric model M of the underlying structural system, estimating the corresponding unknown parameters θ. However, the use of either excessively simplified or too complex models may have a detrimental effect on the possibility to track the system state: oversimplified models may underestimate the effect of a physical process taking place; on the other hand, complex models may lead to good data fitting but possibly yield to poor predictions. In the latter case, the model overfits the incoming data. In [2], an online model class selection strategy was proposed in the framework of an EKF's parameter estimates. Here, a similar approach was adopted for a simultaneous parametric estimate and model class selection exploiting the UKF. Adopting a number n m of possible model classes, the model evidence (or plausibility) consisting of the probability p(M m k+1 ) ∈ (0, 1) was computed for each model class M m , with m = 1, . . . , n m at each time step t k+1 . The sum of the n m model evidences is equal to the unity. To derive the expression of p(M m k+1 ), first, the Bayes theorem was used, giving where p y k+1 |M m k , called conditional evidence, represents the contribution of the measurement at t k+1 to the plausibility of the m-th model class.
Second, we extended the procedure explained in [2] from the EKF to the UKF. As a result, at the end of the corrector phase (after Step 15 of Algorithm 1), the following expression for the conditional evidence applied: where det(·) calculates the determinant of the input matrix. The reported expression approximates p y k+1 |θ, M m k due to the use of Laplace's asymptotic expansion [2].  rhs k+1 q k ,q k ,q k ,θ k+1|k+1 displacement field 15: Computeq k+1 andq k+1 velocity and acceleration fields 16: end for

Results and Discussion
As a numerical case study, we studied how to determine the interstorey stiffness and damping of the two DOF shear building models (n = 2) reported in Figure 1. The mechanical properties of the building were adimensionalised to ease the UKF tuning by setting the matrices in Equation (1) equal to The building was excited by the ground acceleration a 0 reported in Figure 2, lasting 60 s. The response of the building was monitored by recording the floor acceleration y = [y 1 , y 2 ] T with a sampling frequency of 50 Hz for a total of n t = 3000 samples. A white noise, featuring a standard deviation of 5 × 10 −3 , was added to y 1 and to y 2 to mimic the signal perturbation affecting micro-electro-mechanical accelerometers [12].  The acceleration recordings coming from this reference building were used as measurements in the corrector phase of the filtering procedure (Steps 9 and 10 of Algorithm 1). Three model classes, M 1 , M 2 and M 3 , featuring different structural parametrizations, were considered, as shown in the following: Model class M 1 is governed by the parameter θ 1 1 ruling the interstorey stiffness of both floors (for this reason, θ 1 1 is factored out from K 1 ); M 2 is governed by θ 2 = θ 2 1 , θ 2 2 T ruling, respectively, the interstorey stiffness and damping of both floors; M 3 is governed by θ 3 = θ 3 1 , θ 3 2 , θ 3 3 T , where θ 3 1 and θ 3 2 rule the first and second floor interstorey stiffness and θ 3 3 rules the damping associated with both floors. Comparing these parametrizations with the reference model, it is clear that M 1 is underparametrizing the mechanical system, not associating any parameter with the damping properties of the structures and suffering a model bias, being C 1 = 0.6 C; M 3 is overparametrizing the stiffness matrix and M 2 is performing a correct parametrization of the structural response, and it is therefore expected to allow for the best estimate of the system mechanical properties. For all model classes, the initial guesses of the relevant parameters underestimated by 40% of the parameter values ruling the reference structure.
KF tuning is usually problem-dependent and is performed through a trial-and-error procedure. In this case, we have set the SP scaling parameters to α SP = 10 −3 , κ SP = 0 and β SP = 2; the measurement noise covariance to R = 4 × 10 −4 I 2 , where I 2 ∈ R 2×2 is the identity matrix; the process noise covariance to Q = 10 −8 I n p , with I n p ∈ R n p ×n p ; and the initial parameter covariance toP θθ 0 = 0.25 I n p . The value of n p depends on the number of parameters employed by each model (n p = 1 for M 1 , n p = 2 for M 2 and n p = 3 for M 3 ).
In Figure 3, the predicted output of M 1 , computed according to Step 8 of Algorithm 1, is reported against the floor acceleration measurements, showing the filter capacity of tracking the shear building accelerations despite the presence of noise. A small discrepancy between the reference model and the predicted output is observable only magnifying the curves. The predicted outputs of M 2 and M 3 , not reported for lack of space, exhibit an even smaller discrepancies. The filter capacity of tracking the system output was expected to greatly help parameter identification. In Figures 4-6, the time evolution of the parameters employed by M 1 , M 2 and M 3 are reported, respectively. Black colour is used for parameters involved in the expression of the structural stiffness; orange colour when related to the structural damping. The plots report both the parameter posterior estimates and the confidence intervals of these estimates. Looking at the confidence intervals, stiffness-related parameters seem to assume negative values during the first part of the analyses. This is due to to the initial choice ofP θθ 0 = 0.25 I n p . However, positive values have been always associated with the interstorey stiffness due to the use of the scaled version of the UKF. Similar reasoning applies to damping-related parameters.  Looking at Figure 4, the UKF was unable to provide a correct estimate for θ 1 1 , despite the uncertainty reduction linked to the narrowing of the confidence interval. Even the stiffness-related parameters θ 3 1 and θ 3 2 of M 3 , depicted in Figure 6, seem not able to converge to the desired value. On the contrary, coming to M 2 , θ 2 1 was correctly identified with small uncertainty, as shown in Figure 5. These results were somehow expected due to the underparametrization of the mechanical system operated by M 1 and the overparametrization of the mechanical system exhibited by M 3 , while M 2 embodied the correct description of the reference model.
Model class M 3 was unable to provide any idea of the damping properties, ending up pushing θ 3 3 to 0. Model class M 2 provided a better estimate, still quite poor, overestimating by 40% the damping related parameter θ 2 2 . These difficulties were due to the relevance of damping in the identification of continuously excited structures, discussed in [4].
From the results reported above, M 2 seems to lead to the best system identification; however, we reached this conclusion by knowing the mechanical properties of the reference system. It would have been very hard, if not impossible, to judge model plausibility simply looking at the predicted outputs. Indeed, as shown in Figure 3, the UKF has been able to reproduce the monitoring system outcome even when M 1 is employed. For this reason, model evidence computation, whose outcome is reported in Figure 7, is extremely relevant to understand which model can be trusted the most. At the beginning of the identification procedure, equal plausibility was associated with the three models. Their values were recursively updated as soon as new measurements became available using Equations (7) and (8). During the first part of the analysis, M 1 appeared to be the most plausible model class. This is in agreement with intuition: M 1 is the easiest to tune, employing just one parameter, and the bias in the modelling of damping has a marginal relevance when t < 20 s due to the strong ground motion undergone by the structure. In a second stage, M 3 resulted to be the most plausible model class. This was due to the good estimate of both the stiffness-related parameters and the dampingrelated parameters in the central part of the analysis. Finally, the overcomplexity of M 3 led to a deterioration of the parameter identification, while the good convergence of the stiffness-related parameters and the reasonable damping estimate promoted M 2 as the most plausible model class.
This numerical example shows that model evidence evaluation can be successfully used for model selection. The reader should note that, due to the recursive nature of Equation (7), a certain time delay occurred between the improved identification capacity of the filter equipped with a certain model and the increase in plausibility of this model.

Conclusions
In this work, we have discussed an algorithm for simultaneous parameter estimation and model evidence calculation in dynamic linear elastic problems. Starting from the work of [2], a recursive expression for model evidence evaluation was derived for when the unscented Kalman filter is used. Numerical results show that model evidence can guide system identification in the presence of model uncertainties by associating a plausibility measure with different employed models featuring possible parametrizations of mechanical domains. Indeed, model evidence can be successfully used to select the most plausible structure parametrization as parameter identification is carried out.