Identification

Practical Aspects of System Identification

2. Selection (Estimation) of the Model Complexity

This section discusses the determination of the delay (d) and the degrees of the polynomials A(q^-1), B(q^-1) and C(q^-1).

The “plant + disturbance” model to be identified is of the form


A(q ^1 ) y(t)  q ^d B(q ^1)u(t)  w^(t)

where, according to the structures chosen, one has:


w^(t)  e(t)


w^(t)  A(q ^1 )w(t)
(7.3.1)

(7.3.2)

w^(t)  C(q ^1)e(t)
(7.3.3)

w^(t) 
1


C(q ^1 )


e(t)

(7.3.4)

In order to start the parameter estimation methods, we need to specify

n_A  ? ;
n_B  ? ;
d  ?

and for structures [S3] and [S4]

n_C  ?

where n_A, n_B and n_C are the degrees of the polynomials A(q^-1), B(q^-1) and C(q^-1) respectively.

Remember that the order of a sampled-data system is given by n = max (n_A, n_B+d).
Techniques for order estimation have been presented in Chapter 6, Section 6.5. If order estimation functions that allows the estimation of n = max (n_A, n_B+d), n_A, n_B and d, are available (like estimorderiv.sci (Scilab) and estimorderiv.m (MATLAB^®)¹, that estimate n), it is advisable to use them. Nevertheless the order estimation can also be done on the basis of “a priori” knowledge and a set of trials on the acquired data². We give next details about this procedure.

“A priori” Choice for n_A

Two cases can be distinguished:

1. 
   Industrial plant (temperature control, flow rate, concentration, etc...).
   

For this type of plant in general


n _A  3
¹ Available from the book website.
² This approach can be also used for a confirmation of the order estimation results obtained with the techniques described in Section 6.5.

and the value n_A= 2, which is very typical, is a good starting value to choose.

2. 
   Electromechanical systems.
   

n_A results from the structural analysis of the system.


Example: flexible transmission with two resonant modes.
In this case, n_A=4 is chosen, since a second-order is required to model each resonant mode.


Initial Choice of d and n_B
If no knowledge of the time delay is available, d = 0 is chosen as an initial value. If a minimum value is known, an initial value d = d_min is chosen.
If the time delay has been underestimated during identification, the first coefficients of B(q^-1) will be very low. Thus n_B must then be chosen so that it can indicate the presence of the time delays and identify the transfer function numerator. n_B = (d_max - d_min) + 2 is then chosen as the initial value. At least two coefficients are required in B(q^-1) because of the fractional delay (see Section 2.3.7). If the time delay is known, n_B  2 is chosen, however the value 2 is a good initial value.


Initial Choice of n_C
As a general rule n_C = n_A is chosen.


Determination of the Time Delay d (First Approximation)
Method 1. One identifies the system using the RLS. The estimated numerator will be of the form


B^ˆ(q ^1 )  b^ˆ q ^1  b^ˆ q ^2  b^ˆ q ^3  ...
1 2 3

If

1
b^ˆ  0.15 b^ˆ₂


b₁  0 is considered, and time delay d is increased by 1: d = d_in + 1 (since if b₁ = 0, B(q^-1) = q^-1 (b₂ q^-1 + b₃ q^-2)).
If

i
b^ˆ_i  0.15 b^ˆ_{d 1}
; i  1,2,...,d_i

time delay d is increased by d_i: d = d_in + d_i. After these modifications, a new identification session is performed.
Method 2. An identification of the system is performed with the RLS in order to obtain a model of the type “impulse response”


n_B
y(t)  b_iu(t i) ;
i1
n_B  Large (20,...,30)

If there is a time delay, then

b

i
^ˆ  0.15 b^ˆ_{d 1}
; i  1,2,..., d

From this time delay estimation, a new identification is performed to find a pole- zero model.

Both these methods may of course be completed by a display of the step response.
A more accurate estimation of time delay d is carried out by performing a new identification, followed by a validation of the identified model. Note that if the system is contaminated by measurement system noise an accurate estimation of the delay will only be carried out with the method that enables the identified model to be validated.


Determination of (n_A)_max and (n_B)_max
The aim is to obtain the simplest identified model that verifies the validation criteria. This is linked on the one hand to the complexity of the controller (which will depend on n_A and n_B), and, on the other hand, to the robustness of the identified model with respect to the operating conditions.
From the results presented in Section 6.5, a first approach to estimate the values of (n_A)_max and (n_B)_max is to use the RLS and to study the evolution of the variance of the residual prediction errors, i.e. the evolution of

R(0)  E ²
(t)
1 _^N

 ² (t)

^N t 1
as a function of the value of n_A + n_B. A typical curve is given in Figure 7.5.
In theory, if the example considered is simulated and noise free, the curve should present a neat elbow followed by a horizontal segment, which indicates that the increase in parameter number does not improve the performances. (see Section 6.5, Figure 6.5.1).
In practice, this elbow is not neat if non-white noise is present. The practical test used for determining n_A + n_B is the following: consider first n_A, n_B and the corresponding variance of the residual errors R(0). Consider now n_A' = n_A + 1, n_B and the corresponding variance of the residual errors R'(0).

If


R^(0)  0.8R(0)
it is unwise to increase the degree of n_A (same test with n_B' = n_B + 1).

R(0)

The good value

1 2 3 4 5 6 7

^{npa=r (nA+ nB})

Figure 7.5. Evolution of the variance of residual errors as a function of the number of model parameters

It is advisable to increase the values of n_A and n_B in a non-simultaneous way. A simultaneous increase may induce the identification of a pair of very close poles and zeros that have a small influence on the criterion with respect to the increase of n_A or n_B but that could make the controller computation harder.

With the choice that results for n_A and n_B, the model identified by the RLS does not necessarily verify the validation criterion. Therefore, while keeping the values of n_A and n_B, other structures and methods must be tried out in order to obtain a “valid” model. If after all the methods have been tried none is able to give a model that satisfies the validation criterion, then n_A and n_B must be increased.
The estimated numerical values of the coefficients corresponding to the highest

powers of polynomials
A^ˆ(q^1) ^and B^ˆ(q^1) ^{(after the estimation of the time delay}

1. 
   also give a clue on the maximum order to be chosen. If these values are very small compared to the previous ones, it is often possible to reduce the order of A^ˆ(q^1) ^and B^ˆ(q^1) ^{(one must obviously perform a new identification session with}
   

the updated orders).
The order estimation techniques presented in Section 6.5 rapidly provide good results. Moreover, it is recommended to compare these results with the available a priori information upon the system structure.

2. Identification of Simulated Models: Examples³

We consider now two files (T and T1) containing I/O data generated by a known discrete-time model. Each file contains 256 input/output recordings. In each case, the input has been a PRBS generated by a register with 8 cells (N = 8). The length of a complete sequence is therefore 255⁴.

The file T has been obtained using the following model:


A(q ^1) y(t)  q ^1B(q ^1)u(t)

where
A(q ^1)  1  1.5q ^1  0.7q ^2
B(q ^1)  1q ^1  0.5q ^2