A New Maximum Likelihood Estimator Formulated in Pole-Residue Modal Model

Recently, a lot of efforts have been devoted to developing more precise Modal Parameter Estimation (MPE) techniques. This is explained by the necessity in civil, mechanical and aerospace engineering of obtaining accurate estimates for the modal parameters of the tested structures, as well as of determining reliable confidence intervals for these estimates. The Non-linear Least Squares (NLS) identification techniques based on Maximum Likelihood (ML) have been increasingly used in modal analysis to improve precision of estimates provided by the Least Squares (LS) based estimators when they are not accurate enough. Apart from providing more accurate estimates, the main advantage of the ML estimators, with regard to their LS counterparts, is that they allow for taking into account not only the measured Frequency Response Functions (FRFs) but also the noise information during the parametric identification process and, therefore, provide the modal parameters estimates together with their uncertainties bounds. In this paper, a new derivation of a Maximum Likelihood Estimator formulated in Pole-residue Modal Model (MLE-PMM) is presented. The proposed formulation is meant to be used in combination with the Least Squares Frequency Domain (LSCF) to improve the precision of the modal parameter estimates and compute their confidence intervals. Aiming at demonstrating the efficiency of the proposed approach, it is applied to two simulated examples in the final part of the paper.


Introduction
Over the last years a lot of efforts have been spent in developing more accurate modal parameter estimation techniques. The outcome of these efforts are the development of precise and robust estimators, as, for instance, the poly-reference Least Squares Complex Frequency Domain (pLSCF) [1,2] formulated in frequency domain and the Stochastic System Identification (SSI) algorithms developed to estimate the modal parameters from the measured responses [3]. Initially formulated for input-output systems and also known by its commercial name PolyMAX [2], the pLSCF was afterwards extended to output-only systems [4]. One of the main advantages of this technique is the possibility of creating clear stabilization diagrams, since the physical poles tend to stabilize faster over the identified model orders when compared to other Modal Parameter Estimation (MPE) algorithms, as, for instance, the poly-reference Least Squares Complex Exponential (pLSCE) (also known as LSCE-Prony) [5]. Later on, the capability of estimating confidence intervals for modal parameter estimates was added both to the LSCF [6] and pLSCF [7] estimators, as well as to the SSI [8][9][10] technique.
Recently, new improvements were added to the estimation with the pLSCF algorithm to overpass the loss of precision and robustness observed when dealing with very noisy Frequency Response Functions (FRF)s and poorly excited modes [11][12][13]. It is verified that, under these circumstances, the algorithm is not so robust and tends to underestimate the damping ratios. This new variant of the algorithm consists of a 2-step approach which was developed to address these issues [14]. In fact, this approach consists of a combination of two different estimators, for example, the pLSCF and the poly-reference Maximum Likelihood Estimator in Modal Model formulation (pMLE-MM) [12,[14][15][16][17] used in a final step of the identification process to increase the precision of the estimates provided by pLSCF as well as to compute the uncertainties bounds on the estimated modal parameters.
The combination of these two estimators enhances the MPE by taking advantage of the specific features of each one of them, as for instance, the fast and clear stabilization diagrams provided by the pLSCF and the precision and numeric stability of the pMLE-MM to improve the accuracy of the identified modal parameters and estimate their confidence intervals. More recently, new improvements were added to the combined pLSCF-pMLE-MM by imposing identification constraints based both on the assumption of reciprocity of the identified modal model and of normal (real) mode shape vectors [18]. In this paper, a new ML Estimator (MLE) formulated in Pole-residue Modal Model (MLE-PMM) is proposed to improve the accuracy of the modal parameters estimated with the Least Squares Complex Frequency (LSCF) [1]. Similarly to the combined pLSCF-pMLE-MM, the idea of using the MLE-PMM in combination with the LSCF estimator is to retain and improve multi-reference information provided by the latter in a Non-linear LS sense (NLS) with the former. Similarly to the non-linear optimization technique described in Reference [19], the proposed ML-based approach is formulated with the modal model in pole-residue form.
The innovative characteristics of the latter with regard to the former, however, are related to the fact that: (1) both the FRFs and the noise information are taken into account in the modal parameter estimation and, therefore, apart from providing the modal parameter estimates the technique also provides the uncertainties on the estimated modal parameters; (2) the initial values can be calculated by the LSCF estimator that gives good starting values in easy way by means of a clear stabilization chart and (3) the optimization is done by using an advanced technique that enlarges the convergence region and prevents the algorithm from converging to local minima. The paper is basically divided in three parts. The LSCF identification technique and the Least Squares Frequency-Domain (LSFD) [20] estimator in pole-residue modal model with upper and lower residual terms are briefly described in the first part. Next, the novel MLE-PMM is derived in the the second part and, finally, in the last part of the paper, the implementation of the MLE-PMM herein proposed is evaluated, both in terms of accuracy of the modal parameter estimates and precision of the predicted confidence bounds, using two simulated examples.

Structure of the Combined LSCF-MLE-PMM
The ML-based technique proposed in this paper optimizes the invariants of the Pole-residue Modal Model iteratively and, therefore, it requires a good starting guess to assure convergence over the performed iterations. Though the model parameter estimates provided by any LS estimator can be used as starting guess, some particular advantages are obtained if the ML-based estimator described in the second part of the paper is combined with the LSCF identification technique formulated in Common Denominator Model (CDM). These advantages include: (1) the fact that the Pole-residue Modal Model and the CDM are closely related and, therefore, the LS multi-reference information provided by the LSCF estimator can be retained and iteratively optimized by the MLE-PMM; and (2) the possibility of combining the ML Estimator formulated in CDM (MLE-CDM) [21] with the LSCF estimator to improve even more the accuracy of the LSCF estimates used as starting guess.
The identification process with the combined LSCF-MLE-PMM consists basically of the following steps: (1) primary identification of the physical poles with LSCF estimator from the measured FRFs and estimation of the corresponding modal residues together with the lower and upper residuals with the Least Squares Frequency Domain (LSFD) estimator [20]; and (2) optimization of the model parameters provided both by the LSCF and LSFD estimators and estimation of their respective uncertainty bounds with the MLE-PMM in the final step of the identification process. The second step of this combined technique is derived in this paper and consists of the novel ML-based estimator itself. The chart flow of the identification process with the combined LSCF-MLE-PMM is shown in Figure 1.

Modal Identification with the LSCF and LSFD Estimators (1st Step)
This step aims at estimating the poles corresponding to the natural frequencies and modal damping ratios of the tested structures with the LSCF technique. This estimation is normally carried out by means of stabilization charts constructed by identifying models with increasing order. Once the poles are identified, they are subsequently used by the LSFD technique to estimate the modal residues in a LS sense. For completeness, a brief description of both techniques used in the first step of the identification process with combined LSCF-MLE-PMM is presented in following subsections.

The LSCF Estimator
The first derivation of the LSCF aimed at identifying input-output systems using the measured FRFs as primary data [5,20]. Afterwards, this derivation was extended to the identification of output-only systems [22]. This technique estimates the invariants of the Common Denominator Model which is assumed to represent the measured FRFs and is given bŷ where N i and N o denote the number of inputs and outputs, respectively, with the numerator matrix N k (Θ, ω) and the denominator scalar d(Θ, ω) given, respectively, by and In Equations (2) and (3) the parameters Ω r (ω) = z r represent discrete-time polynomial base, with z = e −jω∆t denoting the z-domain variable, ∆t the sampling time and ω the circular frequency. The coefficients α r and β kr in Equations (2) and (3) are the unknown invariants to be estimated. Detailed discussion about the LSCF technique can be found in several publications and, therefore, it is not discussed in the context of the present paper. Further details about the implementation of this estimator can be found, for instance, in References [5,20,22].

The LSFD Estimator
The LSFD is a technique commonly used in frequency-domain modal analysis to estimate the modal residues and to model the influence of the out-of-band modes in the frequency band of interest. When the poles are known a priori, these parameters can be estimated in a LS sense with the following pole-residue modal model with upper and lower residual terms where H(ω) ∈ C N o ×N i is the modelled FRF matrix; N o , N i and N m denote, respectively, the number of outputs, the number of inputs and the number of vibration modes; λ m ∈ C and [Res] m ∈ C N o ×N i stand for the pole and modal residue corresponding to the mth vibration mode; LR, UR ∈ C N o ×N i are, respectively, the lower and upper residuals used to model the influence of the out of band modes in the considered frequency band and a and b are scalars defined according to the measured quantity (i.e., displacement, velocity or acceleration) and to the type of analysis performed (i.e., Experimental or Operational Modal Analysis (EMA or EMA)) in Table 1 [4]; and the operator (•) * denotes the complex conjugate and j = √ −1 the imaginary unit. The poles λ m occur in complex-conjugated pairs and are related to the eigenfrequencies ω n m and damping ratios ξ n m as where ω n m = 2π f n m , with f n m representing the natural frequencies of each mode. An extended description of the LSFD estimator formulated in Pole-residue Modal Model is found in Appendix A. Once the poles are estimated with the LSCF technique and the modal residues with the LSFD algorithm, these parameters can be optimized by means of the ML estimator formulated using the Modal Model in Pole-residue form. The idea behind the MLE-PMM algorithm is to optimize the initial modal parameters given, for instance, by the LSCF and the LSFD estimators by adjusting these values so that FRF equation in Pole-residue Modal Model (4) is fitted to the measured FRF. Compared to the LS-based approaches with uncertainty estimation capabilities, the ML-based algorithms have the advantage of taking into account both the measured FRFs and the noise information at once during the parametric identification process and thus, apart from providing the more accurate estimates for the modal parameters, they also allows for computation of confidence intervals for these estimates directly from the Fisher information matrix.
Furthermore, the ML approach herein introduced is formulated in continuous-time with the Pole-residue Modal Model, meaning that it does not suffer from numerical ill-conditioning which is a well-known identification problem of the rational fraction polynomial models in Laplace domain. The optimization of the starting parameters with the MLE-PMM is accomplished by minimizing the following (negative) log-like cost function with E k (Θ, ω f ) denoting the equation error given as follows: where H k (ω f ) and σ H k (ω f ) are, respectively, the measured FRF and its standard deviation and H k (Θ, ω f ) is the modelled FRF; k (k = 1, 2, · · · , N o N i ) denotes any of the N o N i elements of the measured and estimated FRF matrices; N f represents the number of frequency lines and ω f = 2π f the angular frequency evaluated at frequency line f . The parameter Θ is a column vector containing the invariants of the pole residue model to be optimized by the ML algorithm and is given as follows with and The ML estimates of the invariants of the modal model is given by minimizing the cost function (6) with respect to parameter Θ. Similarly to the MLE formulated in Commom Denominator Model (MLE-CDM) described in Reference [21], this minimization is accomplished by means of a Gauss-Newton optimization algorithm combined with Levenberg-Marquardt approach to guarantee a continuous reduction of the cost function over the performed iterations. In case of real coefficients, each Gauss-Newton iteration is performed in two steps:

1.
Solve the normal equations

2.
Compute an update of the previous solution the Jacobian matrix and ∆Θ i the perturbation on the parameters Θ evaluated at the ith iteration. The error equation or so-called residual vector is given by . . .
and the Jacobian matrix by ; and X k ∈ C N f ×2N m the derivatives with respect to the natural frequencies and damping ratios. The entries Y k and X k are computed, respectively, by (15) and with the partial derivatives in each entry computed according the expressions found in Appendix C. By taking advantage of the block structure of the Jacobian matrix, the normal Equation (11) are rewritten as follows with

Fast Estimation of the Perturbations on Modal Parameters
From Equation (17), the perturbation on the coefficients ∆θ k (i.e., the perturbations on the real and imaginary parts of the modal residues and on the real and imaginary parts of the upper and lower residuals) can be written as a function of the perturbation on the poles, ∆θ λ , as By making use of Equation (18), the perturbations ∆θ k can eliminated from the normal Equation (17), yielding or in a more compact form This elimination decreases the memory, as well as the time required to run the algorithm. As for the MLE-CDM [21], an efficient implementation of the MLE-PMM is only possible if the variances are taken into account in the cost function (6). Also, because it uses the noise information as a weighting in its cost function, the MLE-PMM can handle measurements with a large dynamical range. Once the perturbations on the natural frequencies and damping ratios are calculated by means of Equation (19), then perturbations on the modal residues and on the upper and lower residuals are computed using Equation (18).
It is worth highlighting that, although the derivation of the reduced normal equations as in Equation (20) follows a strategy similar to those in References [1,14], the resulting system of equations is different. This follows from the fact that the system of equations represented by Equation (20) is formulated with the Pole-residue Modal Model with upper and lower residuals (4). The derivation of the reduced normal equations as in Equation (21) is described in Appendix B.

Estimation of the Covariance of the Measured FRFs
The estimation of the uncertainty bounds on the identified modal parameters is only possible if the variance of the noise is taken into account in the cost function (6). In the formulation of the MLE-PMM proposed in this paper, it was considered the following assumptions for the estimation of this variance [7]: • the noise on the measured FRFs is circular complex normally distributed; • the noise on the measured FRFs is zero-mean valued; • the noise on the measured FRFs is uncorrelated over the frequencies; and • the noise on the measured FRFs is uncorrelated over the outputs.
The estimation of this variance depends on the estimator used in the non-parametric estimation of the FRFs. If the H 1 and H 2 estimators are used to estimate the FRFs, the variance can be computed from the coherence function [22]. In case of Single Input Multiple Output (SIMO) systems, the variance of the H 1 FRF estimator can be calculated by [23] (22) where N b is number of time data segments used to estimate the FRF and γ 2 oi is the coherence function. This function is commonly used in EMA to indicate how the output is correlated with the input. The closer to one this function is, the more the outputs are influenced by the inputs. On the other hand, the closer to zero it is, the more the outputs are contaminated by noise. Assuming Multiple Input Single Output (MISO) and Multiple Input Multiple Output (MIMO) systems with noise uncorrelated over outputs, the covariance matrix of the oth row of the FRF matrix is given by [22] Cov where γ 2 o (ω f ),Ŝ y o y o andŜ uu are the multiple coherence function, the autopower spectrum of the output and the autopower spectrum of the input, respectively, with the subscript y o denoting the DFT spectrum of a single measured output.

Convergence of the ML Algorithm
The MLE-PMM algorithm optimizes iteratively the starting estimates provided by the LSCF. This is accomplished by setting a maximum number of iterations and a value for the relative error between two consecutive iterations of the ML-PMM algorithm (i.e., |l(Θ i+1 ) − l(Θ i )| /l(Θ i ), with l denoting the cost function to be minimized and i the iteration number). In order to avoid local minima and enlarge the convergence region of the Gauss-Newton algorithm, the following Levenberg-Marquardt form of Equation (11) should be used [23] to ensure that the cost function decreases during the performed iterations Increasing the Levenberg-Marquardt parameter λ LM in Equation (24) forces the cost function to decrease but reduces the convergence speed. Normally λ LM = 0 is used as an initial value and then it is adapted in every iteration according to the variation of the cost function. If the cost function increases, the value of parameter λ LM is increased, otherwise, it is decreased [22].

Estimation of the Uncertainty Bounds
One of the main advantages of the ML-based algorithms is the possibility of estimating the confidence intervals for the identified modal parameters using the noise information measured together with the FRFs during the vibration tests. As shown in Reference [24], a good approximation of the covariance of the ML parameters Θ MLE−PMM is obtained by with J l standing for the Jacobian matrix evaluated in the last iteration of the Gaussian-Newton algorithm. By taking advantage of the structure of the Jacobian matrix and using the matrix inversion lemma [25], the covariance of the natural frequencies and damping ratios can be estimated independently from the covariance of the modal residues and from the covariance of the lower and upper residuals, as follows The advantage of the parametrization used in Equation (8) is that, since all the linearization is implicitly performed during the computation of the Jacobian matrix, the variances of the natural frequencies and damping ratios are computed directly from the reduced normal matrix, M 1 . Therefore, the use of additional linearization formulas is avoided if one is only interested on these variances. However, if one is also interested on the variances of the real and imaginary parts of the poles, they need to be linearised around Cov( f n , ξ n ) by means of the following linearization formulas The covariance of the modal residues and the upper and lower residuals are estimated by Once the covariances of the real and imaginary parts of the modal residues are computed, the covariance of the mode shapes and of the modal participation factors are estimated by following the procedure described in Reference [26].

Logarithmic MLE-PMM
Another variant of the proposed MLE-PMM can be derived by minimizing the following "log-like" cost function with σ log H k (ω f ) denoting the standard deviation of the kth element of logarithmic FRF matrix. The relationship between this standard deviation and the measured one is given by which is a measure of the noise-to-signal (power) ratio of a complex logarithmic function [27].
The strategy used to minimize the cost function (6) can be employed to Equation (29) and the partial derivatives of the logarithmic error Equation (30) with respect to the invariants of the pole-residue modal model are computed with the expressions presented in Appendix D.

Simulated Data Analysis
In order to validate the proposed implementation of the MLE-MM, two examples are analysed. The first consists of a simulated EMA of a five-DOF system and the second corresponds to a lattice tower structure with two pairs of closely spaced modes. In these analyses, the modal parameters together with their uncertainty bounds were estimated using the new MLE-PMM introduced in Section 2.2. The results obtained for both examples are presented in following subsections.

Five-DOF System
The first example used to validate the ML-based approach discussed in Section 2.2 is illustrated in Figure 2. This system was used in Reference [28] to compare different modal parameter estimation techniques in terms of their sensitivity to statistical errors. It is composed by 5 masses supported by cantilever beams and connected among themselves by arch springs. The exact natural frequencies, damping ratios and modal masses of the system are given in Table 2, whereas the real modes are shown in Table 3. These properties were used to generate the FRFs used in simulated EMA. The system was excited by a white Gaussian noise at masses 1 and 2 and the responses were measured at all DOFs, resulting in FRF matrix with two columns and five rows. These FRFs were calculated in the frequency range of 0-80 Hz with a resolution of 0.1 Hz. Afterwards, a colored noise was introduced in the FRF matrix with a standard deviation of 10%. The noise was added to the real and imaginary parts independently and was calculated as a percentage of the absolute value of the FRF at each frequency line. This was achieved by adding a complex random number to the FRF at each frequency line. This number was computed so that its amplitude is a random number of a normal distribution and its phase is an uniform random number between 0 and 2π [13].  The exact and noisy element(1,1) of the FRF matrix and the corresponding "exact" standard deviation of the noise are shown in Figure 3a. A set of 5000 FRFs with 800 frequency lines contaminated with noise was generated to perform Monte Carlo simulations in order to assess the efficiency of the proposed MLE-PMM. The modal parameters of each dataset were identified with the LSCF and LSFD estimators and then used as starting values to be optimized by the MLE-PMM algorithm. The identification of each dataset was performed using the full frequency band, that is, with no upper and lower residual terms. Though 1000 Monte Carlo simulations would be more than enough to assess the robustness of the MLE-PMM, a higher number of realizations were performed, in this particular simulation example, to also verify its efficiency in predicting the 95% elliptical confidence regions. A typical stabilization diagram constructed with the LSCF method from the noisy FRF is shown in Figure 3b.    Analyzing these figures, one verifies that, while no changes can be visually perceived for the estimated natural frequencies with regard to their exact counterparts, a significant improvement in accuracy of the damping ratio estimates after the optimization with MLE-PMM is clearly visible. This conclusion is corroborated by the results presented in Table 4 where the results obtained from the Monte Carlo simulations are compared in terms of mean relative bias error. One of the main advantages of the MLE-PMM approach herein proposed with regard to the optimization technique described in Reference [19] is that the former takes into account not only the FRFs but also the noise information as primary data and, therefore, apart from the modal parameter estimates, it also predicts the uncertainties on theses estimates. In order to illustrate this capability, the standard deviations of the natural frequencies and damping ratios predicted after 10 iterations of the proposed MLE-PMM for the 3rd and 5th modes in each Monte Carlo run are compared to the respective sample standard deviations in Figures 6 and 7. In these figures, the results obtained with the proposed MLE-PMM are also compared to those from the MLE-MM described in Reference [14]. Similarly to the MLE-PMM, the optimization of the pLSCF estimates with the MLE-MM was also performed with 10 iterations in each Monte Carlo run. It is clear from Figures 6 and 7 that the standard deviations provided by the proposed MLE-PMM are in very good agreement both with their corresponding sample standard deviations and with standard deviation obtained with MLE-MM.
The estimates for the damping ratios obtained after 10 iterations of the MLE-PMM in each Monte Carlo run are plotted as a function of the corresponding natural frequency estimates in Figure 8. Also, in this figure, the 'classical' 95% confidence ellipses calculated both with the sample covariance and the mean covariance computed from the estimates obtained with Equation (26) over 5000 Monte Carlo runs are shown. It is clearly seen from such figure that not only the shape of the uncertainty regions defined by the damping ratio and natural frequency estimates but also the estimated 95% confidence ellipses are in good agreement with the corresponding sample covariance ellipses.    These results also demonstrate the MLE-PMM capability of providing very accurate confidence regions for the optimized modal properties. The MLE-PMM capability of predicting the sample standard deviations for all the identified modes is illustrated in Figure 9 and Table 5. In the former, the sample standard deviations computed from the 5000 sets of estimates obtained with the Monte Carlo simulations are compared to the means of the standard deviation estimates provided by the MLE-PMM. In the latter, this comparison is established, in terms of relative difference, both for the MLE-PMM and the MLE-MM estimates. Analysing these results, one verifies that the MLE-PMM provided accurate predictions for the sample standard deviation of all the identified modes, which are as well in good agreement with the MLE-MM estimates.  Although it is verified from Figures 6 and 7 and from Table 5 that the MLE-MM provides slightly lower uncertainty bounds for optimized estimates, the MLE-PMM algorithm is significantly less time consuming. This is verified in Figure 10 where the performance of the ML-PMM algorithm is compared, in terms of time consumed over the performed Monte Carlo simulations, to its ML-MM counterpart. The significant difference in processing time between the two algorithms is explained by the fact that the dimensions of the system matrices obtained with the MLE-PMM is significantly lower than those of the MLE-MM.

Lattice Tower Example
The second simulated example used to validate the ML approach discussed in Section 2.2 consists of a finite element model which is illustrated in Figure 11. The example consists of a lattice tower structure constituted by two segments of equal height and by variable equilateral triangular sections. The lower section is scaled with regard to the upper one by a factor of 2. This model is composed by beam elements with 6 DOFs per node. The nodes of the foundations are clamped and the others have 3 DOFs: two translations in x and y and one rotation around z-axis; and the remaining DOFs are set equal to zero. These settings result in a FE model with a total of N m = 18 DOFs.
The symmetry of the structure was slightly broken by defining different geometry and materials to the column elements in order to simulate the behaviour of tower-like structures, which normally present very closely spaced of bending modes. The real mode shapes of the tower structure are shown in Table 6. The first 6 modal configurations and the corresponding modal parameters of the FE model of the lattice tower are presented in Figure 12 and Table 7, respectively.
The structural damping, on the other hand, is modelled as the special case of proportional damping by setting the damping coefficients of all modes equal to 1%. The structure was excited, independently, at node 8 in x and y-direction and at node 9 in y-direction with white noise inputs. The simulated responses were measured in acceleration at the 6 DOFs indicated in Figure 13. By assuming that the two triangular sections behave as rigid diaphragms, the measured outputs are enough to yield the modal configurations of the structure. The FRF matrix with 6 rows and 3 columns was calculated in the frequency range of 0-10 Hz with a resolution of 10 mHz. Afterwards, a colored noise was introduced in the FRF matrix with a standard deviation of 5%.    The noise was added to the FRFs following the same procedure described in the previous example. The exact and noisy element(1,1) of the FRF matrix, the corresponding exact standard deviation of the noise as well as the exact natural frequencies are shown in Figure 14. A set of 1000 FRFs with 1000 frequency lines contaminated with noise was generated to perform Monte Carlo simulations in order to assess the efficiency of the proposed MLE-PMM when dealing with closely spaced modes. In the first step of the identification of each dataset, a stabilisation diagram was constructed with the LSCF technique from the measured FRFs and the approach based on Hierarchical Clustering described in References [29,30] was subsequently used to automatically identify the physical modal properties.
A typical stabilization diagram constructed by identifying models with order ranging from 30 to 60 is represented in Figure 15. Next, the FRF matrix was synthesized from the estimates obtained with the LSCF technique was then used as starting guess by the MLE-PMM. Finally, a total of 20 iterations were performed, both to optimize the estimates obtained with LSCF in the first step and to compute confidence intervals of the optimized estimates. Differently from the previous application example, a higher number of iterations was used in these simulations to secure the convergence of the cost function l(Θ). This difference is explained by the fact that the system being optimized in this particular simulation example comprises vibration modes with closely spaced natural frequencies and, therefore, the ML estimator needs to iterate more times to reach convergence. Similarly to the previous application example, the optimization with the MLE-PMM was carried out using again the full frequency band, i.e., with no out-of-band modes.    The variation of the natural frequency and damping ratio estimates for modes 2 and 5 obtained both with the LSCF and MLE-PMM over the Monte Carlo simulations are shown in Figures 16 and 17 and the results obtained in terms of relative bias error for all the identified modes are synthesized in Table 8.  By comparing the results obtained before and after the optimization with MLE-PMM shown in this table, it is verified that no significant reduction in terms of the bias error was obtained after 20 iterations of the MLE-PMM with regard to the LSCF estimates. The mean of standard deviation estimates for the natural frequencies and damping ratios of modes 2 and 5 predicted by the proposed MLE-PMM after 20 iterations are compared to their respective sample standard deviations in Figures 18 and 19. By inspecting these figures, one verifies that, even in case of closely spaced modes, the MLE-PMM provides very good estimates for the sample standard deviations. These results are also verified in Figure 20 and Table 9, where the sample standard deviations of all the identified natural frequencies and damping ratios are compared to those predicted by the MLE-PMM. In Figures 18 and 19 and in Table 9, the estimates for the standard deviations provided by the MLE-PMM are also compared to the estimates obtained with MLE-MM proposed by Reference [14] after 20 iterations. From this comparison it is concluded that, even though the MLE-PMM provides very accurate estimates for the uncertainty bounds of its optimized modal parameter estimates, these uncertainty intervals are significantly higher than those provided by the MLE-MM. This leads to the conclusion that, in case of closely spaced modes, the optimisation with the MLE-PMM is not as efficient as in case of reasonably or well separated modes.

Conclusions
In this paper, a novel ML-based estimator formulated in pole-residue modal model is proposed. Just like any other ML-based approach, the MLE-PMM introduced in this paper requires a good starting guess to assure convergence over the performed Gauss-Newton iterations. Despite the fact that any LS-based identification technique can be used to provide the starting guess, the proposed ML approach is combined with the LSCF estimator to retain the multi-reference information and improve the accuracy of the LS estimates. The idea of combining the ML-PMM with LSCF estimator aims at taking advantage of the particular features of each estimators, that is, (i) the fast modal parameter estimation obtained with the latter; and (ii) the numerical stability of the former. It was verified that, when combined with the LSCF estimator, the MLE-PMM is capable of improving the accuracy of the LSCF estimates and providing confidence intervals for the optimized estimates.
The main advantage of the proposed MLE-PMM with regards to the MLE-MM is that, due to the reduced dimension of the system matrices, it is significantly less time consuming, as demonstrated in the first application example presented in Section 3.1. Another distinguishing feature of the derived ML approach is that, due to the parametrization used in its formulation, the variances of the natural frequencies and damping coefficients are estimated directly from the Jacobian matrix. Therefore, the use of additional linearization formulas is avoided if one is only interested on these variances. The main disadvantages, however, are in fact, related to the limitations of modal model in pole residue form. These disadvantages include: (1) difficulty of distinguishing between close spaced modes; (2) difficulty in accurately synthesizing the FRF matrix once the mode shapes and modal participation factor vectors are estimated from the modal residues by means of the singular value decomposition; and (3) difficulty in fully optimizing the LSCF modal parameter estimates used as starting guess in case of closely spaced modes.
When dealing with reasonably or well separated modes, however, it is observed that MLE-PMM performs very well when compared to it modal model based counterpart. Two application examples were used to validate the proposed ML-based estimator. The first consisted of well separated whereas the second of closely spaced modes. From the simulations with former, it was demonstrated the ability of the MLE-PMM technique to significantly reduce the bias error of the LSCF estimates and predict the uncertainties for the optimized estimates. This is demonstrated by the results presented in Tables 4 and 5. From the simulations with the latter, it was verified (as shown in Tables 8 and 9) that, even in case of closely spaced modes, the proposed MLE-PMM was capable of reducing the bias error of the LSCF estimates and providing reliable predictions for the confidence intervals of the optimized estimates.
with the subindex "o" designating the oth row of a matrix, E(Θ, ω) the error matrix and Θ representing the row vector containing the unknown parameters given by Writing down Equation (A1) for each frequency line f ( f = 2, 3, · · · , N f ), yields with E o (Θ) and H o given now by and the corresponding Jacobian matrix, J ∈ C (N f −1)×(2N m +2) , by It is straightforward to solve Equation (A3) for Θ in a linear least squares sense. Yet, in order to guarantee the realness of Θ, J and H o are replaced in Equation (A3), respectively, by and, after some manipulations of the resulting expression, the following equation can be derived to compute the modal residues and the upper and lower matrices of the residual model

Appendix B. Derivation of the Reduced Normal Equations
Since the Jacobian matrix as in Equation (14) normally contains a huge amount of elements, the computation of the inverse of the normal equations as in (17) slows down the estimation process with the MLE-PMM. Furthermore, the numerical stability of the Gauss-Newton algorithm can be compromised depending on the conditioning of this matrix. In order to avoid these issues, the reduced normal equations can be computed instead. The derivation of the reduced normal equations for the MLE-PMM described in Section 2.2 starts by rewriting Equation (17) in a more compact form, as From Equation (A8), one can write Isolating ∆θ in the first set of Equations (A10) and inserting the resulting expression into the second set, yields Now replacing Equations (A9) into (A11), gives

Appendix C. Partial Derivatives of the Classical ML-PMM
The partial derivatives in entries Y k and X k of the classical ML-PMM, given by Equations (15) and (16), are computed by . . .

Appendix D. Partial Derivatives of the Log-Like ML-PMM
The partial derivatives in entries Y k and X k of the log-like ML-PMM, given by Equations (15) and (16), are computed by . . .