The Stochastic Subspace Identification Techniques

 

By

 

Palle Andersen                                                    Rune Brincker

Structural Vibration Solutions A/S                 Aalborg University

NOVI Research Park                                           Department of Building Technology and Structural Engineering

Niels Jernes Vej 10,                                            Sohngaardsholmsvej 57,

DK-9220 Aalborg East.                                     DK-9000 Aalborg

Denmark                                                                Denmark

 

1 Introduction

 

In the traditional input-output modal analysis the estimation of modal parameters have been performed using a somewhat deterministic mathematical framework. One of the major hurdles for people of this traditional modal community to overcome, when turning to output-only modal analysis, is the switch of the mathematically framework. In output-only modal analysis the mathematically framework involves the use of statistics and introduction of concepts such as optimal prediction, linear system theory and stochastic processes.

 

The two general assumptions made in output-only modal analysis is that the underlying physical system behaves linearly and time-invariant. The linearity imply that if an input with a certain amplitude generates an output with a certain amplitude, then an input with twice the amplitude will generate an output with twice the amplitude as well. The time-invariance implies that the underlying physical system does not change in time. One of the typical parametric model structures to use in output-only modal analysis of linear and time-invariant physical systems is the stochastic state space system.

                                                                                                                                  (1)

The first part of this model structure is called the state equation and models the dynamic behavior of the physical system. The second equation is called the observation or output equation, since this equation controls which part of the dynamic system that can be observed in the output of the model. In this model of the physical system, the measured system response yt is generated by two stochastic processes wt and vt. These are called the process noise and the measurement noise. The process noise is the input that drives the system dynamics whereas the measurement noise is the direct disturbance of the system response.

 

The philosophy is that the dynamics of the physical system is modeled by the n´n state matrix A. Given an n´1 input vector wt, this matrix transforms the state of the system, described by the n´1 state vector xt, to a new state xt+1. The dimension n of the state vector xt is called the state space dimension. The observable part of the system dynamics is extracted from the state vector by forward multiplication of the p´n observation matrix C. The p´1 system response vector yt is a mixture of the observable part of the state and some noise modeled by the measurement noise vt.

 


2 The Statistical Framework

2.1 Properties of stochastic state space systems

 

The state space model (1) is only applicable for linear systems that do not have time-varying changes of its characteristics. However, this is not the only restriction. The only way to obtain an optimal estimate of a state space model on the basis of measured system response, is to require that the system response is a realization of a Gaussian distributed stochastic process that has zero mean.

 

In other words, in the applied stochastic framework the system response is modelled by a stochastic process yt defined as

                                                                                                                                (2)

and the principal assumption is that the measured system response is a realization of this process. It is seen that this process is completely described by its covariance function. This means that if we can estimate a state space model having the correct covariance function this model will completely describe the statistically properties of the system response. An estimated model fulfilling this is called covariance equivalent. The estimator that can produce such model is called an optimal estimator.

 

Since the system response of the linear state space model is a Gaussian stochastic process it implies that ,  and all are Gaussian stochastic processes as well. Since the input processes  and are unknown we make the simplest possible assumption about their statistical properties, which is to assume that they are two correlated zero-mean Gaussian white noise processes, defined by their covariance matrices as

                                                                                                   (3)

The Gaussian stochastic process describing the stateis also zero-mean and completely described by its covariance function

                                                                                                                              (4)

Using (1) to (4) the following relations can be established

                                                                                                              (5)

The matrix G is the covariance between system responseand the updated state vector. The covariance function of can also be expressed in terms of the system matrices as

                                                                                              (6)

There are two additional system matrices turns out to play an important role

                                                                                                  (7)

These are the extended observability matrixand the reversed extended stochastic controllability matrix .

2.2 Optimal prediction

 

One of the most important parts of all estimation is the ability to predict the measurements optimally. In output only modal analysis this means to be able to predict the measured system response optimally. An optimal predictor is defined as a predictor that results in a minimum error between the predicted and measured system response. If the system response can be predicted optimally it implies that a model can be estimated in an optimal sense.

 

Recall that the state vector  completely describes the system dynamics at time t. In order to predict the system response  optimally it is necessary to start by defining an optimal predictor of . Now assume that we have measurements available from some initial time k = 0 to k = t-1. Collect these measurements in a vector

                                                                                                                                     (8)

In the Gaussian case the optimal predictor of  is then given by the conditional mean value

                                                                                                                                 (9)

So, the optimal predictor of  is defined as the mean value of  given all measured system responsefrom k = 0 to k = t-1. The difference between  and  is called the state prediction error  and is defined as

                                                                                                                                     (10)

This error is the part of  that cannot be predicted by .

 

In order to predict the system response a similar conditional mean can be formulated for

                                                                                                                  (11)

The last part of this equation is obtained by inserting (1) and assuming that andfrom k = 0 to k = t-1 are uncorrelated.

2.3 The Kalman filter.

 

The two predictors (9) and (11) are related through the so-called Kalman filter for linear and time-invariant systems, see e.g. Goodwin et al. [6]

                                                                                                                             (12)

The matrix  is called the non-steady state Kalman gain and is called the innovation and is a zero-mean Gaussian white noise process. Defining the non-steady-state covariance matrix of the predicted state vectorasthe Kalman gainis calculated from

                                                      (13)

The last of these equations is called the Ricatti equation. The Kalman filter predicts the state on the basis of the previous predicted state and the measurement. The covariance                                                                      of the innovations can be determined from the last equation in (12) as

                                                                                                                      (14)

Given that the initial state prediction isand the initial state prediction covariance matrix  and assume that we have measurements available from k = 0 to k = t-1, then this filter is an optimal predictor for the state space system (1) when the measurementsare Gaussian distributed.

2.4 The innovation state space system.

 

At start up the Kalman filter (12) will experience a transient phase where the prediction of the state will be non-steady. However, if the state matrix A is stable the filter will enter a steady state as time approach infinity. When this steady state is reached the covariance matrix of the predicted state vectorbecomes constant, i.e. , which imply that the Kalman gain becomes constant as well, i.e.. The Kalman filter is now operating in steady state and is defined as

                                                                                                                               (15)

The steady state Kalman gain is now calculated from

                                                              (16)

The last equation is now called an algebraic Ricatti equation. Assuming all matrices but P is known this equation can be solved using eigenvalue decomposition, see Aoki [2] and Overschee et al [1].

 

If the last equation in (15) is rearranged the following state space system is obtained

                                                                                                                            (17)

This system is called the innovation state space system. The major difference between this system and (1) is that the state vector has been substituted with its prediction, and that the two input processes of (1) have been converted into one input process – the innovations. This state space system is widely used as model structure in output only modal analysis, see e.g. Ljung [3] and Söderström et al. [4].

 

 

3 The Stochastic Subspace Identification Framework

 

The Kalman filter defined in the last section turns out to be the key element in the group of estimation techniques known as the stochastic subspace techniques.

 

From (17) it is seen that if sufficiently many states of (1), let’s say j states, can be predicted, i.e.  and , then the A and C matrices can be estimated from the following least regression problem

                                                                           (18)

This is a valid approach since the innovations are assumed to be Gaussian white noise. Since A and C are assumed to be time-invariant this regression approach will be valid even though the predicted state  and  originates from a non-steady state Kalman filter.

 

So the fundamental problem to solve in the stochastic subspace identification technique is to extract the predicted states from the measured data. To show how this is performed, consider the state space system in (1) and take the conditional expectation on both sides of both equations to yield

                                                                                         (19)

Now assume that a recursion is started at time step q. Inserting the first equation in (19) recursively into itself i times and finally inserting the result the last of the equations in (19) leads to the following formulation

                                                                                                                                    (20)

This equation shows the relation between the initial predicted state and the prediction of the free (noise free) response of the system . By stacking i equations on top of each other the following set of equations are obtained

                                                                                             (21)

By introducing the vector as the left-hand side and noticing that the first part of the right-hand side is equal to the extended observability matrix we actually obtain the following expression for the predicted states

                                                                                                                                        (22)

The matrixis actually the pseudo-inverse of . This equation shows that if we can estimate and for several values of q, we can in fact estimate the predicted states for several values of q as well.

3.1 Estimation of free system response.

 

In this section we will focus on the estimation of the predicted free response . We will estimate a set of vectors and gather them column by column in a matrix O.

 

In order to predict the system response a conditional mean similar to (11) can be formulated.

                                                                                                                      (23)

This conditional mean is the prediction of the future system response  given the past system response from time t = i+q-1 down to t = q. This conditional expectation is only an approximation of (11) since the conditioning vector stops a time t = q and not t = 0. The approximation is only good if i is sufficiently high. For zero-mean Gaussian stochastic processes this conditional mean can be calculated by, see e.g. Melsa et al. [5].

                                                                       (24)

Since the error is zero-mean and uncorrelated and is independent of the conditional mean and the conditioning vector the conditional mean (24) is also called the orthogonal projection of the vector onto the vector .

 

In order to estimate all elements of we need to extend (24) to allow estimation of  to in one operation. This is done by using (8) to extend the conditional mean  in (24) to the following . This results in the following equation for

                                                           (25)

In the last equation a new ip´ip matrix  is introduced for simplicity. This matrix is defined as

                                                                                                    (26)

Incidentally, the matrix  is also equal to

                                                                                                                                          (27)

As seen in (18) we need a bank of predicted state estimates for q = i to q = i+j-1 for a sufficiently large value of j. To estimate these state in one operation based on the approach in (25) we need to define the following two matrices and as

 

                                                                       (28)

 

                                                                           (29)

 

The index p in (29) signifies that the matrix contains system response of the past compared to the system response we are predicting. Since we assume that the system response is stationary, i.e. that