Mechanical Systems and Signal Processing (1999) 13(5), 803}819 Article No. mssp.1999.1224, available online at http://www.idealibrary.com on
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS FOR MODAL ANALYSIS M. SMAIL AND M. THOMAS Ecole de Technologie SupeH rieure, Montreal, Quebec, Canada AND
A. A. LAKIS Ecole Polytechnique, Montreal, Quebec, Canada (Received 25 June 1998, accepted 23 February 1999) The Autoregressive moving average (ARMA) model is a very e$cient technique for modal parameter identi"cation of mechanical systems, especially when the signal is noisy. However, when signi"cant noise is present in the signal, it is necessary to increase the computational order of the ARMA model. Consequently, this arti"cial increase of the model order yields to a more di$cult identi"cation of the exact number of modal parameters in a given frequency range, especially when we have no prior knowledge of the behaviour of the mechanical system. A new method based on the eigenvalues of a modi"ed covariance matrix is proposed. It is shown that the eigenvalues of the covariance matrix that lead to a minimum and constant value depending on the noise level, correspond to supplementary orders induced by the noise. Thus, the exact order of the mechanical system is revealed from the analysis of the eigenvalue magnitudes with the model order. The analysis of the gradient of the eigenvalue computed at the exact order allows also to select the minimal and necessary order used for computation, without any prior modal parameter identi"cation. This method is robust to noise level and sensitive to the sampling frequency. Thus, the application of the proposed method at di!erent sampling frequencies allows to select the optimal sampling frequency by reducing the lack of accuracy in the identi"cation of modal parameters. 1999 Academic Press
1. INTRODUCTION
The Autoregressive moving average (ARMA) or ARX way of modelling has generated great interest in the broad "eld of identi"cation and evaluation of parameters in such di!erent areas as signal processing, speech analysis, image processing, telecommunications, spectral analysis, automation, etc. There are three important stages in the analytical process: selection of the model, identi"cation and evaluation. A wide variety of methods and techniques for solving the di!erent problems have been put forward and developed. Once a model has been selected on the basis of observations or prior knowledge of the signal, the use of these methods requires basic knowledge of the systems to be analysed, especially of the model order to be used and the evaluation of the parameters obtained. The selection of the computing model order is the "rst step in the identi"cation process. Several authors [1}4] have discussed the question and have proposed certain criteria. Among the better known are the criteria based on the statistical properties of residual prediction errors, the "nal prediction error based on the variance of the prediction error, which was proposed by Akaike (AIC) [5], and the minimum description length (MDL) developed by Schwartz and Rissanen [6]. The application of these criteria requires "rst the 0888}3270/99/050803#17 $30.00/0
1999 Academic Press
804
M. SMAIL E¹ A¸.
selection of a possible interval for the model order to be used, for a series of pairs (p, q), and then an evaluation of the parameters. The selection of the order is made on the basis of minimal variance. It is necessary, therefore, to have a prior evaluation with di!erent orders, but this results in signi"cant computation time. More recently, other authors [7, 8] have proposed new approaches derived from the MDL criterion and based on the minimal singular values of the covariance matrix without prior evaluation. To determine the model order, it is necessary to determine a limiting upper value in order to establish the order of the matrix. However, it is di$cult to determine that value when signi"cant experimental noise is present. More speci"cally, in the mechanics of vibration, certain criteria are used, such as the stability diagram which involves noting the parameters which stabilize as the computing model order rises. The modal assurance criterion (MAC) [9] requires the calculation of the mode shapes and the consideration of these with coherence close to unity. The modal assurance criterion provides a measure of consistency between estimates of a modal vector. It can only indicate consistency, not validity. It is a scalar constant relating the causal relationship between two modal vectors, taking a value from zero to one. Hence, if the modal vectors under consideration exhibit a consistent relationship, the MAC should approach unity. The MAC can be used to verify or correlate an experimental modal vector with respect to a theoretical modal eigenvector. It could be used also to evaluate modal parameter estimation methods if a set of analytical frequency response functions with realistic levels of random and bias errors is generated and used in common to a variety of modal parameter estimation methods. For using these two previous methods, a prior evaluation is necessary. It is obvious that the application of these criteria requires a lengthy period of calculation and investigation. Knowledge of the exact order of the model is necessary but not su$cient. In the evaluation of the parameters with an exact order model, it is not possible to obtain all the modes if the measurements are corrupted by experimental noise [10 }12]. To obtain all the modes with acceptable accuracy, a higher order is required [13, 14]. Therefore, two orders have to be established before beginning the identi"cation procedure: the exact order which corresponds to the number of real frequencies and the necessary order which will allow us to obtain acceptable accuracy when noise is present. In this paper, we present a new method for establishing the exact model order and the necessary model order for acceptable accuracy. A great advantage of this method lies in the fact that no prior evaluation of the parameters is necessary. The method is applied to the AR and ARMA models.
2. THE ARMA PROBLEM
The general formulation of an ARMA or ARX model [15] appears as follows: A(z\)y(k)"B(z\)u(k) where y(k) is de"ned as the output signal and u(k) as the input signal with A(z\)"1#a z\#a z\#2#a z\N N B(z\)"b z\#b z\#2#b z\O O and z\ is the operator such that z\y(k)"y(k!1).
(1)
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
805
Whether the series u(k) is visible or not, it is called ARX or ARMA. The stability condition requires that the solutions of the polynomials A(z\) lie within the unit circle. Furthermore, it is taken that the input u satis"es the condition of sustained excitation. The identi"cation in this case is based on the evaluation of parameters a and b which are G G the unknown quantities of the problem. For a selected order p and q, the model can be written in the following way: y(k)"!a y(k!1)!a y(k!2)!2!a y(k!p)#b u(k!1) N #b u(k!2)#2#b u(k!q) O which can be rewritten in the form of the product of vectors: y(k)"+u(k!1),2+h,
(2)
(3)
where +u, is the observation vector: + (k!1),2"+!y(k!1),!y(k!2), 2 ,!y(k!p), u(k!1), u(k!2), 2 , u(k!q), and +h, is the parameter vector +h,"+a , a , 2 , a , b , b , 2 , b ,2. N O This relation is veri"ed at di!erent times, and can be expressed by the following system: y(k)"+u(k!1),2+h, y(k#1)"+u(k),2+h, y(k#2)"+u(k#1),2+h, $ y(k#N)"+u(k#N!1),2+h,
(4)
which can be written in the vector and matrix form: +>,"[Q] +h,
(5)
where +>, is a N;1 vector, +h, is a (p#q);1 vector and [Q] is a N;( p#q) matrix. At this level, we can obtain an estimate by the least-squares method [16]: +hM ,"[Q2Q]\[Q]2+>,.
(6)
It is well known that methods of identi"cation based on the least-squares technique lead to biased estimates once the disturbance is no longer white noise. However, the way in which the ARMA model for an identi"cation problem is formulated leads to the inclusion of a disturbance which is not a white noise in equation (6). Thus, given y(k) the measured value, we can write y(k)"yJ (k)#w(k)
(7)
where yJ (k) is the value without noise and w(k) is white noise with zero mean and variance p, and u(k) is supposed not to be corrupted by measurement noise. Replacing this in the equation for the ARMA model, we have y(k)!w(k)"!a (y(k!1)!w(k!1))!a (y(k!2)!w(k!2))!2!a (y(k!p) N !w(k!p))#b u(k!1)#b u(k!2)#2#b u(k!q). (8) O
806
M. SMAIL E¹ A¸.
By reorganizing the terms expressing the actual measurements and the experimental noise, we obtain y(k)"!a y(k!1)!a y(k!2)!2!a y(k!p)#b u(k!1)#b u(k!2)#2 N #b u(k!q)#m(k) O
(9)
with m(k)"w(k)#a w(k!1)#a w(k!2)#2#a w(k!p) (10) N where the disturbance is no longer white noise and consequently the evaluation will be biased.
3. IDENTIFICATION METHOD USING CORRECTED COVARIANCE MATRIX METHOD
We return to the ARMA model as expressed by equation (3): y(k)"+u(k!1),2+h,.
(11)
We split the part of the observation vector formed by the output signal from that formed by the input signal, and we can then write +u(k!1),2"(+u (k!1),2; +u (k!1),2). (12) W S In this vector, the assumption is made that only y(k) is corrupted by measurement noise (output error method). In the same way, we separate the parameters of the AR part of the parameter vector from those of the MA part: +h,2"+a , a , 2 , a , b , b , 2 , b ,"++a,, +b,,2 N O
(13)
with +a,2"+a , a , 2 , a , N +b,2"+b , b , 2 , b ,. N Taking these notations into account, the ARMA model [17] can be expressed as follows: y(k)"+uy(k!1),2+a,#+u (k!1),2+b,#+m(k),. (14) S In the following, we examine the e!ect of experimental noise on the covariance matrix. We start by considering the relationship between the estimate obtained by the least-squares method and the covariance matrix starting from equation (6): [C]+hM ,"+X,
(15)
where [C] is a matrix ( p#q, p#q) and +X, is a ( p#q);1 vector: +u , +u ,2 ! +u , +u ,2 1 C !C W W W S WW WS [C]"[Q]2[Q]" " I I N ! +u , +u ,2 !C C +u , +u ,2 SW SS S W S S I I
(16)
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
807
and
or more simply as
1 y(k) +uW, +X,"[Q]2 +>," I N y(k) +u , S I
(17)
1 [C]" +u(k), +u(k),2 N I 1 [X]" y(k) +u(k),. N I
(18)
As we have assumed that the excitation is not corrupted by the measurement noise and also the w(k) is white noise with zero mean and variance v*, we can show here the e!ect of the measurement noise on the covariance matrix. The matrix components [C ], [C ] and [C ] of [C] denote the classical Toeplitz WW WS SS covariance matrix of signals y and u, where [C ] can be written as follows: WW r (0) r (1) r (p!1) WW WW WW [C ]" r (!1) r (0) r (p!2) WW WW WW WW r (1!p) r (2!p) r (0) WW WW WW
(19)
where 1 1 1 r (0)" y(k)" (yJ (k)#w(k))" y(k)#v* WW N N N I I 1 1 r (i!j)" y(k!i);y(k!j)" (yJ (k!i);yJ (k!j)) WW N N I and the same for [C ], [C ] and +X,. WS SS Since the measurements include noise and it is presumed that this is white noise, it is easy to see that matrix [C] can be written as
[C ]"E[C]"
C J J #v*I !C J WW WS !C J C SW SS
(X )"E(X)
(20) (21)
whereyJ (k)"y(k)!w(k), the value not contaminated by experimental noise and E de"nes the expectation. If we de"ne [C ] as a function of v, [C ] having been calculated as above using the T T contaminated values y(k),
C !vI !C WS . [C ]" WW T !C C SW SS
(22)
The parameter vector +h , is de"ned by T [C ] +h ,"+X, T T
(23)
808
M. SMAIL E¹ A¸.
with
a +h ," T . T b T The identi"cation problem of the ARMA model becomes one of "nding the &&best v'' to satisfy equation (15). We de"ne a mean square error as follows, for large N: 1 D (v)" (+u,2+h ,!y(k)). T , N I By clarifying parts AR, MA and y(k)"yJ (k)#w(k) in this equation, D(v)"E(+u ,2+a,#+u , +b,!yJ (k)!w(k)) W S we have the ARMA model which satis"es the following equation:
(24)
(25)
yJ (k)"+u J ,2+a,#+u ,2+b, S W
(26)
D*"E(=!w(k))"( #a ##I)v* T
(27)
from which
with ="a w(k!1)#a w(k!2)#2#a w(k!p). N By developing equation (24), D(v)"+h ,2[C ] +h ,!2+h ,2 +X,#p T T T W
(28)
where +h ,2 [C ] +h ,"+h ,2[C ] +h ,##a #v. T T T T T T Thus, including equation (29) in equation (28), we obtain
(29)
D(v)"!+h ,2 +X,#p##a #v. (30) T W T It is clear that the &&best v'' is obtained for D(v)"D*; thus equations (27) and (30) give d(v)"D(v)!D*"!+h ,2 +X,#p!v (31) T W which is called an &&error'' function. The expression d(v) will be replaced by its approximate form for large N leading to the equation 1 d(v)"!+h ,2 +X,# y(k)!v. (32) T N Note that d(v) depends explicitly on v and indirectly, non-linearly, on v through the intermediary h . T The study of the asymptotes of d(v), as indicated earlier, leads to a solution for all the null eigenvalues of C . This can be shown, following equation (23), by T C !v I !C a 0 WW ? WS ? " . (33) !C C b 0 SW SS ?
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
809
By eliminating b , equation (33) becomes ? [([C ]![C ] [C ]\[C ])!v I] +a ,"+0, WW WS SS SW ? ?
(34)
which appears now as an eigenvalue problem. The general appearance of the error function is shown in Fig. 1, where the v (i"1, 2 , n) are the asymptotes which equal the singular values calculated in equation ?G (34). The asymptotic values "x the boundaries for the search of solutions which will cancel out the error function corresponding to the optimal values of the experimental noise variance. In the case of an AR model, the problem (34) leads to the calculation of the eigenvalues of the matrix C . From a study of the eigenvalues of this matrix, the rank of the matrix can be WW established and hence the associated model order. The low singular values correspond to the supplementary orders induced by the noise. On the other hand, the calculation of the singular values, in the case of an ARMA or ARX model, is based on the matrix de"ned by equation (34), in which we "nd the AR part of the matrix C , the term for the MA part C and also the cross-terms AR and MA C . By WW SS WS considering all these previous terms used in equation (34) that allow to take into account the in#uence of the MA terms and the cross-terms AR and MA, the calculation of the singular values would give important information on the real order of the model. By considering the AR or ARMA model, the e!ective presence of modes in a signal will be expressed by singular values which are relatively more signi"cant compared to those due to experimental noise. 3.1. CASE OF MULTIVARIATE MODELS In the case of multidimensional signals (SIMO), the procedure can be applied in two di!erent ways. One way to approach the problem is to perform the identi"cation globally
Figure 1. Variation of the function d(v).
810
M. SMAIL E¹ A¸.
for the whole set of ARMA parameters, using a common noise variance v* for all channels. In this case, all signals correlation contribute simultaneously to the global identi"cation. Alternatively, the identi"cation can also be applied separately on each channel with independent noise variance for each signal. This second approach will result in di!erent estimates of AR parameters and an averaging procedure must be de"ned to select the AR parameters. This procedure will help to select existing modes for modal analysis. Assuming di!erent noise variances on di!erent channels is more realistic from an experimental point of view because each signal is emitted from a di!erent sensor with its own perturbation.
4. ASSESSMENT OF THE ARMA MODEL ORDERS AND GRAPHIC REPRESENTATION
First, the eigenvalues of matrix (34) are computed and ordered in decreasing order. Normalisation with respect to the highest value is carried out in the following way: v k" ?G v ? for i"1, 2, 2 , p. The evolution of the ratio k is represented as a function of di!erent orders, in which each curve corresponds to a given order. A three-mode example of the evolution of k with the model order is shown in Fig. 2.
4.1. SELECTION OF THE REAL SYSTEM ORDER Once the real order of the system is found, the value of k becomes stable to a minimum value no matter the order used. If the order is not satisfactory, the k value does not reach a stable level. To select the exact model order, we notice the order at which all the curves lead to a constant value. This indicates the values from which the other eigenvalues do not
Figure 2. 3-dof system (noise: 10%; sampling frequency: 145 Hz):
, 6; , 8;
, 10;
, 12;
, 14;
, 16.
811
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
vary or vary only slightly. In general, this level begins with the "rst uneven eigenvalue after the exact order. For example, Fig. 2 describes a 3-dof case where the values became stable after the model order of 7, so it is easy to deduce that the real model order would be 6 even if the real order were not known. 4.2. SELECTION OF THE ARMA COMPUTING ORDER The selection of the minimum necessary ARMA computing order is based on the shape of the curve. Of course, the di!erent curves should meet at the common convergence point. The best model order to use for computation is the one which shows the most signi"cant curve slope. The slope is calculated between the common convergence point and the previous one. For example, Fig. 2 shows that the slope computed between the model order of 6 and 7, is more signi"cant for a computing model order of 16.
5. SIMULATION RESULTS: EFFECTS OF NOISE LEVEL AND SAMPLING FREQUENCY ON ACCURACY
In the following subsections, we shall present the results obtained from numerically simulated data and consider the in#uence of di!erent factors. First, we consider the in#uence of experimental noise on singular values of a three-mode system. Then, for a "xed noise level, we analyse the in#uence of sampling frequency on the assessment of the exact and necessary system order, which is a very important aspect of the identi"cation process. For numerical simulation, we consider frequencies and damping e!ects which correspond to the natural frequencies of a beam over 0}100 Hz, which was analysed experimentally elsewhere. The values used are given in Table 1. The "gures which follow represent the singular values normalised by the largest value in each of the di!erent orders selected. 5.1. EFFECT OF MEASUREMENT NOISE The sampling frequency is set at a working minimum (145 Hz) for identifying the third mode and the noise level is varied from 10 to 400%. Figure 2 shows the normalised singular values plotted with the model order for a low noise level (10%) by selecting in the computing ARMA order, orders ranging from 6 (exact order) to 20. We can see that for the di!erent orders selected, the ratio k converges toward a value close to zero after the sixth order, after which the level becomes unchanged for the higher orders. This analysis reveals that the exact order of the system can thus be assumed to be 6. Elsewhere, we have noticed that a higher level of accuracy on modal parameters was obtained for the curve showing the largest gradient of k, computed at the convergence order. This is shown in the "gures by a greater slope. Order 16 gives the greatest accuracy and higher orders does not improve signi"cantly this one.
TABLE 1 Modal parameters of a simulated 3-dof system Modes
Mode 1
Mode 2
Mode 3
Frequencies (Hz) Viscous damping rates (%)
3.81 0.2
16.1 0.12
64 0.1
812
M. SMAIL E¹ A¸.
The same remarks can be made for the cases with higher noise levels. At 25, 50 and 100% noise levels, a similar variation is found and the conclusions still hold. The only di!erence is found for the level of the convergence value of k. We note that, as the experimental noise level rises, the convergence value of k also rise, as shown in Fig. 3 for a noise level of 100%. The e!ect of the noise level can be expressed as an increase of the convergence level. With very high noise levels (200}400%), the exact order of the mechanical system can always be identi"ed close to 6 (order where the singular value becomes stable) while the selection of the computing necessary order (slope of the singular value) is more ambiguous (Fig. 4).
Figure 3. 3-dof system (noise: 100%; sampling frequency: 145 Hz):
, 6;
, 8;
, 10;
Figure 4. 3-dof system (noise: 200%; sampling frequency: 145 Hz): , 20; , 22.
, 8;
, 10; , 12;
, 12;
, 14;
, 14;
, 16;
, 16.
, 18;
813
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
5.2. EFFECT OF SAMPLING FREQUENCY 5.2.1. Model with three degrees of freedom The selection of the sampling frequency is a very important stage in the process of identi"cation and the quality of the evaluation depends upon it. In addition to Shannon's condition, oversampling must be avoided in order to keep a good accuracy in modal parameter identi"cation and an optimal sampling frequency must be found. Unfortunately, selecting an optimal sampling frequency is not an easy task when the investigated frequency is not known before identi"cation is made. To illustrate this point, a noise level of 10% is set and the sampling frequency is varied from 145 Hz (Fig. 1) to 1000 Hz. Figure 5 shows the singular value plotted with the model order for di!erent computing orders and a sampling frequency of 300 Hz. As in Fig. 1, Fig. 5 suggests to select a computing order of 16 in order to "nd an exact order close to 6. However, when increasing the sampling frequency to 500 Hz, it was shown (Fig. 6) that, even though the exact order can easily be established as the convergence is complete after the order 6, the computing order must be increased to 24 in order to obtain a good accuracy. By selecting a sampling frequency which is too high (1000 Hz in this case), a lack of accuracy for identifying the modal parameters was observed, even with a very high order (34) for the ARMA model. In this case, it was not possible to identify either the exact order of the model or the optimal computing order. 5.2.2. Model with xve degrees of freedom To con"rm the previous results, a simulation has been conducted on a 5-dof system and a noise level of 10% by considering that the number of degrees of freedom was not known and by selecting the sampling frequency without any a priori knowledge of its value. The modal parameters used for this example are given in Table 2. For frequencies lower than 200 Hz that did not respect Shannon's condition for all the natural frequencies of the system, the results show clearly that the system order is lower
Figure 5. 3-dof system (noise: 10%; sampling frequency: 300 Hz):
, 6; , 8;
, 10;
, 12;
, 14;
, 16.
814
M. SMAIL E¹ A¸.
Figure 6. 3-dof system (noise: 10%; sampling frequency: 500 Hz): , 24; , 26.
, 12;
, 14;
, 16;
, 18;
, 20;
, 22;
TABLE 2 Modal parameters of a simulated 5-dof system Modes
Mode 1
Mode 2
Mode 3
Mode 4
Mode 5
Frequencies (Hz) Viscous damping rates (%)
16 1
25 1
54 1
63 1
100 1
than 10, but the exact order seems overestimated with respect to Shannon's condition. Figure 7 shows an example of the estimation of the exact order by selecting a sampling frequency of 100 Hz. Figure 8 shows the results for a sampling frequency of 250 Hz that respects Shannon's condition. In this case, it is easy to select a computing order of 24 in order to identify an order of 10. As the sampling frequency is increased up to 500 Hz, the values of system and computing orders remain unchanged. Consequently, it is possible to conclude that the two orders found are the orders to be used to obtain all the modes with good accuracy. This fact is con"rmed by the identi"cation results, compared to the parameters used in the simulation. However, for a sampling frequency of 1000 Hz, we "nd again that oversampling occurred, which leads to a lack of accuracy in the identi"cation of modal parameters and presents a problem in the determination of the order.
6. EXPERIMENTAL INVESTIGATION
The proposed procedure has been con"rmed experimentally by applying the technique to industrial applications. Di!erent excitations have been used such as white noise and impact
815
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
Figure 7. 5-dof system (noise: 10%; sampling frequency: 100 Hz): , 22; , 24; , 26; , 28.
Figure 8. 5-dof system (noise: 10%; sampling frequency: 250 Hz): , 22; , 24; , 26; , 28.
, 10;
, 12;
, 14;
, 16;
, 18;
, 20;
, 10; , 12; , 14;
, 16;
, 18;
, 20;
hammer. Two experimental applications have been considered. The "rst application considers the modal analysis of a beam, where the beam order in a given frequency range is easy to check, but the noise level is unknown. The second application considers the modal analysis of an underground concrete drain pipe where both the pipe order and the noise level are unknown. 6.1. BEAM SUBJECTED TO WHITE NOISE For the active control of vibration in #exible structures [18], an ARMA identi"cation technique has been applied experimentally to a large structure. The test has been carried out
816
M. SMAIL E¹ A¸.
on a rectangular 5;0.5;0.1 m steel beam, inclined at 303 in a hinged-free con"guration (Fig. 9). The beam was excited with a shaker that generates a white-noise excitation. The response of the structure measured by accelerometers located at equal distances from the free end. A 100 Hz low-pass "lter has been applied to measurements. The presence of the "lter for both the response and the excitation will reduce the number of modes to investigate. Figure 10 shows the results for a sampling frequency of 300 Hz. This analysis allows to select a compacting order of 16 in order to identify a system order of 6. In fact an overestimation of the exact order to 8 for sampling frequencies lower than 200 Hz and a lack of accuracy for sampling frequencies greater than 1000 Hz was found. This shows the importance of the selection of the sampling frequency in the identi"cation
Figure 9. Experimental beam con"guration.
Figure 10. Experimental application to a beam (sampling frequency: 300 Hz): , 16; , 18; , 20; , 22; , 24.
, 6;
, 8;
, 10;
, 12;
, 14;
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
817
Figure 11. Response spectrum of a drain pipe (sampling frequency: 800 Hz, N"1024).
process, which is another point to take into consideration. By the way, a modal analysis simulation of the beam has been conducted and the results are shown in Table 1. We can also note the satisfactory similarity between the results obtained by simulation and those obtained by experiment. 6.2. DRAIN PIPE EXCITED BY AN IMPACT HAMMER In order to verify if modal analysis was able to detect cracks in underground concrete drain pipes [19], a classical modal analysis has been conducted in the frequency domain by exciting the pipe with an impact hammer. The signal was too noisy to clearly identify the number of natural frequencies in a given frequency range from this analysis. Figure 11 shows a typical frequency response function recorded in the "eld. From Fig. 11, it can be observed that if the main frequencies are easily identi"ed close to 125, 175 and 225 Hz, it is more di$cult to distinguish real natural frequencies from noise frequencies at lower frequencies. Also, the analysis has been conducted by applying the ARMA technique, although this method needs a sustained excitation. Figure 12 shows the results for a sampling frequency of 800 Hz. The analysis allows to select a computing order of 40 in order to identify a system order of 20. Knowing that the mechanical system has 10 natural frequencies in the 0}400 Hz frequency bandwidth, it is now easier to identify them from the frequency analysis.
7. CONCLUSION
A new approach to the determination of the ARMA model orders is proposed. The technique allows to select a computing minimum necessary order and to identify the exact order of the mechanical system in a given frequency range. This technique can be applied to the AR, ARMA and ARX models and is supported by a graphic representation of the development of the singular values. The simulation examples show that very good results of order were obtained in determining a model by this technique even in the presence of a signi"cant level of measurement
818
M. SMAIL E¹ A¸.
Figure 12. Experimental application to a drain pipe (sampling frequency: 800 Hz): , 14; , 16; , 18; , 20; , 22; , 24; , 26; , 28; , 30; , 32; , 40.
, 8; , 34;
, 10; , 36;
, 12; , 38;
noise. The ARMA method has been found to be very sensitive to the selection of the sampling frequency. Oversampling led to a lack of accuracy and undersampling to an overestimate of the system order. We have observed, with all the di!erent examples, that the shape of di!erent curves allows to appreciate the quality of the selection of the sampling frequency. We have also noted that, once the sampling frequency was no longer appropriate, the appearance of the di!erent curves tended to become confused, a fact which led us to select the order necessary for more precise identi"cation. The experimental applications have proven the e$ciency of the method. An experimental investigation of a beam has shown results similar to those obtained by simulation and the analysis of the pipe allowed to identify the number of natural frequencies while with a frequency domain analysis, it was very di$cult to "nd the di!erent natural frequencies. Otherwise, the principal advantage of this approach lies in the fact that one has no need to make a prior identi"cation but simply a calculation of the singular values. The generalisation of this method to the multivariate ARMA models is not straightfroward, and will require more development. The approach developed in this paper could be used in two ways; model based on input and output signals or on the output signal only. For the "rst case with input and output signals, the singular values of the modi"ed covariance matrix are calculated using equation (34), this shows the contribution of input signal to determine the model order. In the second case, the model contains only the output signal and the terms including the input signal are equal to zero. Hence, the singular values are calculated using the same relation (34) corresponding to the singular values of the classical covariance matrix.
ASSESSMENT OF OPTIMAL ARMA MODEL ORDERS
819
ACKNOWLEDGEMENT
The support of NSERC (Natural Sciences and Engineering Research Council of Canada) through strategic grants is gratefully acknowledged. The authors would like to thank CRIQ (Quebec Industrial Research Center) and Silentec Ltd., for their collaboration in data collection.
REFERENCES 1. M. H. PERROT and R. J. COHEN 1996 IEEE ¹ransactions on Biomedical Engineering 43, An e$cient approach to ARMA modeling of biological systems with multiple inputs and delays. 2. M. NOURI, N. MIKI and N. NAGAI 1994 IEEE ¹ransactions on Biomedical Engineering 43, ARMA model order estimation based on the SVD of the data matrix. 3. E. J. HANNAN 1980 ¹he Annals of Statistics 8, 1071}1981. The estimation of the order of an ARMA process. 4. P. M. DJURIED and S. M. KAY 1992 IEEE ¹ransactions on Signal Processing 40, 2829}2833. Order selection of autoregressive models. 5. M. KAY 1988 Modern Spectral Estimation: ¹heory and Application, Prentice-Hall Signal Processing Series. Englewood Cli!s, NJ: Prentice-Hall. 6. J. RISSANEN 1978Automatica 14, 465}471. Modeling by shortest data description. 7. A. AL-SMADI and D. M. WILKES 1996 Circuits and Systems Connecting IEEE ¹ransactions 505}508. On estimating ARMA model orders. 8. GANG LIANG, D. M. WILKES and J. A. CADZOW 1993 ¹ransactions on Signal Processing 41, 3003}3009. ARMA model order estimation based on the eigenvalues of covariance matrix. 9. EWINS 1986 Modal ¹esting, ¹heory and Practice. Letchworth, England: Research Studies Press. 10. GONTIER, M. SMAIL and P. E. GAUTIER 1993 Mechanical Systems and Signal Processing 7, 45}56. A time domain method for identi"cation of dynamic parameters of structures. 11. R. BEN MRAD and S. D. FASSOIS 1991 Journal of