CN107844627A - It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method - Google Patents
It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method Download PDFInfo
- Publication number
- CN107844627A CN107844627A CN201710875350.6A CN201710875350A CN107844627A CN 107844627 A CN107844627 A CN 107844627A CN 201710875350 A CN201710875350 A CN 201710875350A CN 107844627 A CN107844627 A CN 107844627A
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- msup
- mfrac
- sigma
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/06—Power analysis or power optimisation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Physics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- General Engineering & Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Medicines Containing Antibodies Or Antigens For Use As Internal Diagnostic Agents (AREA)
- Software Systems (AREA)
- Computer Hardware Design (AREA)
- Operations Research (AREA)
- Probability & Statistics with Applications (AREA)
- Geometry (AREA)
- Algebra (AREA)
- Life Sciences & Earth Sciences (AREA)
- Databases & Information Systems (AREA)
- Evolutionary Computation (AREA)
- Other Investigation Or Analysis Of Materials By Electrical Means (AREA)
- Evolutionary Biology (AREA)
Abstract
A kind of only output Time variable structure modal parameter Bayesian Estimation method disclosed by the invention, belongs to Structural Dynamics technical field.AR factor projections into the function space represented by basic function sequence, are drawn functional sequence multi -components time correlation autoregression model FS VTAR by the present invention；FS VTAR parameter is estimated using bayes method, obtains the Posterior distrbutionp of the prior distribution of the response at all moment, parameter matrix and covariance matrix, likelihood function and priori probability density distribution function, calculating parameter matrix and covariance matrix；The model order and functional space exponent number of functional sequence multi -components time correlation autoregression model are determined according to BIC；Time variable structure modal parameter is asked using time-frozen method, the credibility interval of Time variable structure modal parameter is obtained using Monte Carlo simulation, the uncertainty of unknown parameter is estimated in the identification of linear time-varying modal parameters, the robustness to crossing estimation can be improved, reduce to model structure change sensitivity, more accurately handle shorter response data.
Description
Technical field
Time variable structure modal parameter Bayesian Estimation method is only exported the present invention relates to a kind of, more particularly to it is a kind of based on more
The Bayesian Estimation method for being used for only output linearity Time variable structure Modal Parameter Identification of component time correlation autoregression model, category
In Structural Dynamics technical field.
Background technology
In daily life and industrial production, because the change of condition of work, the aging of structure, normal wear etc. are different former
Cause, the modal parameter of many engineering structures can all change with the time, show time varying characteristic.Such as under driving vehicle excitation
Vehicle-bridge system, the fuel mass carrier rocket, the aircraft in flight course under Additional pneumatic effect that change over time, Yi Jike
Variable space flight mechanism of geometry of expansion etc. shows time-varying characteristics.Needs and scientific and technical development with life production,
Gradually require to be designed time-varying engineering structure, malfunction monitoring, vibration control etc., this just needs to have clearly Time variable structure
Understanding, grasps its changing rule, develops the analysis method to Time variable structure.
Therefore, the important method and approach as Time variable structure dynamical property analysis, Time variable structure Modal Parameter Identification
Study one of emphasis as following Structural Dynamics research.Time variable structure Modal Parameter Identification can recognize Time variable structure
Modal frequency, Mode Shape and modal damping, these parameters have important physical significance, can be that the structure of Time variable structure is set
The application of meter, monitoring structural health conditions, structure failure diagnosis, structural vibration control etc. provides strong support.
According to the difference of analysis domain, two classes are broadly divided into currently for the Modal Parameters Identification of Time variable structure：
The first kind is time domain approach, is broadly divided into two classes：Method based on state-space model and based on arma modeling
Method.
In the method based on state-space model, Liu by it is a series of by output response and additional input data group
Into Hankel matrixes carry out singular value decomposition, obtain the modal parameter of linear time varying system.Liu and Deng reduces state sky
Between method to the susceptibility of noise, and apply the method in the experiment of mobile cantilever beam.
In the method based on arma modeling, Petsounis and Fassois propose a kind of arma modeling of time-varying, are used in combination
In the modeling of non-stationary random vibration.Poulimenos and Fassois investigations and comparative studies is a variety of based on TARMA models
The modeling method of non-stationary random vibration, including unstructured Parameters Evolution, random parameter develop and deterministic parameter develops.
Spiridonakos and Fassois proposes a kind of adaptive sequence of function TARMA methods, and TARMA is solved using B-spline function
The coefficient of model, and in the modeling of nonstationary vibration.Yang etc. proposes a kind of vector based on mobile Kriging type functions
TARMA models, it can preferably handle mutation problems.
Second class is time-frequency domain method, is broadly divided into two classes：Imparametrization method and parametric method.
Imparametrization method directly uses time-domain analysis, for example, wavelet transformation, smoothed pseudo wigner ville disstribution,
Hilbert-Huang conversion etc., the parameterized model independent of any system.Ghanem and Francesco proposes that one kind is based on
The discrimination method of small echo, modal parameter is picked out by solving wavelet expansion equation.Roshan-Ghias etc. uses Smoothing Pseudo
Wigner-Ville distribution method, the natural frequency and damped coefficient of identification system.Xu proposes that one kind is based on response signal Gabor
The Modal Parameters Identification of expansion.
Parametric method characterizes system using the parameterized model in system time frequency domain and picks out the parameter of parameterized model,
Such as time-varying centimetre case mold.Zhou etc. proposes that a kind of frequency of two-step least squares estimation method identification Time variable structure and mode are shaken
Type.Zhou etc. is based further on the modal parameter of time-frequency domain maximal possibility estimation identification Time variable structure, and this method is to model order
It is insensitive, and frequency bandwidth can be selected.Louarroudi etc. is based on noise inputs output supervision and proposed for Periodic Time-Varying Systems
A kind of Modal Parameters Identification.Ertveldt introduces time-varying system state simulation of frequency region discrimination method the analysis of aeroelastic flutter
In prediction.
In Time variable structure Modal Parameters Identification is parameterized, the model parameter of parameterized model is at present mostly using most
A young waiter in a wineshop or an inn multiplies and maximum likelihood method for parameter estimation.Least-squares parameter estimation method assumes that modeling error meets Gaussian Profile, and
The problem of can only handling under Gaussian Profile is assumed.Maximum likelihood method for parameter estimation can contemplate different error distribution forms,
The result and the result of least square obtained under Gaussian Profile hypothesis is consistent.Least square and maximum likelihood parameter Estimation
Method belongs to the category of Frequency school, by researcher and business software extensive use.Least square method exists with maximum likelihood method
Poor robustness during treated estimation, and the change to model structure is very sensitive.In addition, existed with least square method and maximum likelihood method
Poor effect when handling shorter response data.There is no a kind of simple and easy method, user is very complicated in application.
In recent years, bayes method has obtained extensive concern due to the characteristics of its is more efficient and more practical.In shellfish
In this method of leaf, generally believe that probability is subjective and can be modified according to effective information.By the unknown ginseng of parameter model
Number regards stochastic variable as, and after in view of other useful informations and data, it is with these parameters that unknown parameter is inferred accordingly
Posterior probability distribution based on.After initial priori is provided, posteriority is may infer that, is carried out more between priori and posteriority
Newly, final Posterior probability distribution is obtained.
The content of the invention
Have for above-mentioned technical problem and bayes method existing for least square method and maximum likelihood method above-mentioned excellent
Gesture, a kind of only output Time variable structure modal parameter Bayesian Estimation method technical problems to be solved disclosed by the invention are to realize
For only output linearity Time variable structure Modal Parameter Identification, and can obtain linear time-varying modal parameters identification in estimate it is unknown
The uncertainty of parameter；Tool has the advantage that：The robustness for crossing estimation can be improved；Reduce to the quick of model structure change
Sensitivity, the shorter response data of more accurate processing.In addition, the present invention is in the case of professional knowledge background is lacked
It can be operated, the modal identification of linear time-varying structure can be widely used in Structural Dynamics engineer applied.
The purpose of the present invention is achieved through the following technical solutions.
A kind of only output Time variable structure modal parameter Bayesian Estimation method disclosed by the invention, first, uses time discrete
Form represent multi -components time correlation autoregression model VTAR, by time-varying AR factor projections to by a series of given basic function sequences
Arrange in the function space represented, draw functional sequence multi -components time correlation autoregression model FS-VTAR expression-form.Its
It is secondary, the model parameter of FS-VTAR models is estimated using bayes method, obtains the response Y at all moment, i.e. joint probability density
Function, the prior distribution of parameter matrix and covariance matrix is obtained, obtain likelihood function and prior probability in bayes method
Density fonction, the Posterior distrbutionp of parameter matrix and covariance matrix is calculated.Again, it is true according to bayesian information criterion
Determine the model order and functional space exponent number of functional sequence multi -components time correlation autoregression model, i.e. naAnd pa.Frozen using the time
Knot method seeks Time variable structure modal parameter, obtains the credibility interval of Time variable structure modal parameter using Monte Carlo simulation, produces
The uncertainty of unknown parameter is estimated in being recognized to linear time-varying modal parameters, and then has and has the advantage that：It can improve
Robustness for crossing estimation；Reduce the susceptibility to model structure change, the shorter response data of more accurate processing.
Also include instructing the structural analysis in Structural Dynamics field with setting using the modal parameters of above step identification
Meter.
According to a kind of described structural modal ginseng for only exporting Time variable structure modal parameter Bayesian Estimation method, obtaining
Number, the vibration frequency range of engineering structure can be obtained, detect whether to meet the requirement such as engineering structure standard and vibration isolation, can instruct
The design of engineering structure.In addition, obtained modal parameters can also be Time variable structure health monitoring, structure failure diagnosis,
The application of structural vibration control etc. provides strong support, is with a wide range of applications and benefit.
A kind of only output Time variable structure modal parameter Bayesian Estimation method disclosed by the invention, comprises the following steps：
Step 1：Linear time-varying structure system is represented using functional sequence multi -components time correlation autoregression model FS-VTAR
System.
Step 1 concrete methods of realizing, comprises the following steps：
Step 1.1：Represented with the form of time discrete shown in multi -components time correlation autoregression model VTAR such as formulas (1)：
In formula,It is NoThe response vector signal of output channel, e [k] are as caused by measurement and model approximation
Unobservable non-stationary residual error, e [k] average is zero, and time-varying covariance matrix isnaFor VTAR model orders
Number, Ai[k] is the i-th rank time-varying AR coefficient matrixes, and k represents k-th of discrete time point.
Step 1.2：By time-varying AR factor projections into a series of function space represented by given basic function sequences.
Shown in given function space such as formula (2)：
In formula, p be function space dimension, fjFor jth rank basis function vector.
By Ai[k] deploys, that is, obtains function spaceThe functional sequence of vectors time-varying autoregressive (FS- of expression
VTAR), as shown in formula (3)：
Shown in x [k] component form such as formula (4)：
In formula, aij,lmFor matrix AijElement, l represent l-th of output channel.
Step 1.3：Draw functional sequence multi -components time correlation autoregression model FS-VTAR expression-form such as formula (5)
It is shown：
Y=Φ θ+ε (5)
In formula,
Wherein, symbol " T " represents transposition computing.
So far, complete to represent linear time-varying structure using functional sequence multi -components time correlation autoregression model FS-VTAR
System.
Step 2：Using the model parameter θ and covariance matrix ∑ of FS-VTAR models in bayes method estimating step 1
Posterior distrbutionp.
Step 2 concrete methods of realizing comprises the following steps：
Step 2.1：Obtain the response Y at all moment, i.e. formula (5) left end item joint probability density function.
Covariance matrix ∑ (k) in formula (1) is definite value, i.e. ∑ (t)=∑.In the viewpoint of Frequency school, unknown ginseng
Number θ and ∑ are definite value, and bayes method thinks that unknown parameter θ and ∑ are stochastic variable.Therefore, responding Y distribution is
Using parameter θ and ∑ as the probability distribution of condition, its probability density function is represented by a multivariate normal distributions, as shown in formula (9)：
Step 2.2：Obtain the prior distribution of parameter matrix θ and covariance matrix ∑.
Parameter matrix θ and covariance matrix ∑ be it is unknowable, it is theoretical based on linear model standard conjugate gradient descent method,
When e [k] in formula (1) is given as normal distribution, under ∑ known case, shown in θ condition prior distribution such as formula (10)：
In formula,It is parameter matrix θ expected matrix,It is covariance matrix,Expression average is a, covariance matrix is b'sTie up normal distribution,Represent Kronecker product.It is right
In positive definite matrix ∑, prior distribution is taken as inverse Wishart distribution, as shown in formula (11)：
In formula,Represent that it is b to obey the free degree0, Scale Matrixes Ψ0The N of distributionoRank positive definite real matrix.
Step 2.3：Obtain the likelihood function and priori probability density distribution function in bayes method.
From formula (9), in the case where θ and ∑ are given, all time of day response Y probability density function, i.e., seemingly
Right function, as shown in formula (12)：
In formula, the mark of tr () representing matrix, exp () represents natural Exponents, and p (x | y) represents the bar of x when known to y
Part probability density function, | | represent determinant.From formula (10) and multivariate normal distributions equation, in the case where ∑ is given
θ priori probability density function p (θ | ∑) is expressed as：
Based on formula (11), the priori probability density function p (∑) of ∑ is expressed as：
Then from formula (13) and formula (14), the joint priori probability density function p (θ, ∑) of θ and ∑ is：
Step 2.4：The parameter matrix θ and the Posterior distrbutionp of covariance matrix ∑ that calculating bayes method obtains.
Because response Y is unrelated with θ and ∑, p (Y | θ, ∑) marginal likelihood is constant, in the case of Y determinations, θ and ∑
Posterior probability density function be expressed as：
In formula,
Contrast (16) and formula (15), obtains θ Posterior distrbutionp：
With the Posterior distrbutionp of ∑：
In formula, parameter b and Ψ are respectively
Step 3：The model order of functional sequence multi -components time correlation autoregression model is determined according to bayesian information criterion
Number and functional space exponent number, i.e. n in formula (4)aAnd pa。
Step 3 concrete methods of realizing is as follows：
According to bayesian information criterion BIC, as model order naWith functional space exponent number paWhen being optimal, BIC value takes
Minimum value.In practical operation, the response signal once obtained is taken, when building the complete functional sequence of vectors as shown in formula (6)
Become autoregression model, and take different naValue and paValue is modeled, to each group of naValue and paValue calculate BIC value, and by its
It is plotted on same figure, chooses BIC and reach corresponding n during minimum valueaValue and paValue, it is selected exponent number.
Step 4：Time variable structure modal parameter and credibility interval are asked using time freezing method.
Step 4.1：The calculating time freezes Time variable structure modal parameter.
Moment point k is freezed using the method for " time freezes ", by VTAR coefficients Ai[k] obtains modal parameter.Using adjoint
Matrix obtains system pole by generalized eigenvalue problem：
(D[k]-λr[k]I)Vr[k]=0 (21)
In formula, λr[k] andBe respectively k moment D [k] r-th of characteristic value and feature to
Amount, D [k] is obtained by AR coefficient matrixes：
" time freezes " the modal frequency f of system is calculatedr[t] and dampingratioζr[t]：
Realize the Modal Parameter Identification for only output linearity Time variable structure.
Step 4.2：The credibility interval of Time variable structure modal parameter is obtained using Monte Carlo method, obtains linear time-varying structure
The uncertainty of unknown parameter is estimated in Modal Parameter Identification, tool has the advantage that：The robustness for crossing estimation can be improved；
Reduce the susceptibility to model structure change, the shorter response data of more accurate processing.
FromSampling, carry out Monte Carlo simulation.Divide for convenience of from probability
ClothSampling, random matrix is constructed by formula (24)Or
In formula, matrix Z each element is satisfied by unit normal distribution, CMAnd C∑It is the cholesky point of M and ∑ respectively
Solution.
∑ is obtained so as to build W by the expectation E (∑) of covariance matrix ∑, according to inverse Wishart distribution, i.e. formula (25)：
By formula (24) and (25), fromObtain the sample that there is same distribution with parameter matrix θ
W, even θ=W.Obtained parameter matrix θ is sampled every time, time frozen structure modal parameter can be obtained according to formula (23), this
Process is designated as single Monte Carlo simulation.Single Monte Carlo simulation step is repeated, carries out Time variable structure modal parameter Bayes
The multiple Monte Carlo simulation of method of estimation, obtain estimating the credibility interval of modal parameter using the method for point estimation afterwards, obtain
The uncertainty of unknown parameter is estimated in being recognized to linear time-varying modal parameters, and then has and has the advantage that：It can improve
Robustness for crossing estimation；Reduce the susceptibility to model structure change, the shorter response data of more accurate processing.When
The average value of varying structure modal parameter is the average value of multiple Monte Carlo simulation result, and variance is multiple Monte Carlo mould
Intend the variance of result.
Also include step 5：The modal parameter that applying step 1 obtains to step 4 identification instructs the knot in Structural Dynamics field
Structure is analyzed and design.
The modal parameter obtained according to step 1 to step 4 identification, the vibration frequency range of engineering structure, detection can be obtained
Whether meet the requirement such as engineering structure standard and vibration isolation, the design of engineering structure can be instructed.In addition, obtained structural modal ginseng
Number can also provide strong branch for the application of the health monitoring, structure failure diagnosis, structural vibration control of Time variable structure etc.
Hold, be with a wide range of applications and benefit.
Beneficial effect：
1st, compared to existing least square method, a kind of only output Time variable structure modal parameter Bayes disclosed by the invention
Method of estimation, bayes method is applied to pick out modal parameter in only output linearity Time variable structure Modal Parameter Identification, adopted
The credibility interval of Time variable structure modal parameter is obtained with Monte Carlo method, obtains estimating in the identification of linear time-varying modal parameters
The uncertainty of unknown parameter can provide the uncertainty of Modal Parameter Identification, and then have and have the advantage that：It can improve pair
In the robustness for crossing estimation；Reduce the susceptibility to model structure change, the shorter response data of more accurate processing.
2nd, a kind of only output Time variable structure modal parameter Bayesian Estimation method disclosed by the invention, can obtain engineering structure
Vibration frequency range, detect whether to meet the requirement such as engineering structure standard and vibration isolation, the design of engineering structure can be instructed.Separately
Outside, the modal parameters obtained can also be the health monitoring of Time variable structure, structure failure diagnosis, structural vibration control etc.
Application strong support is provided, be with a wide range of applications and benefit.
Brief description of the drawings
Fig. 1 is a kind of flow chart of only output linearity Time variable structure Modal Parameters Identification of the present invention；
Fig. 2 is Three Degree Of Freedom spring-dampers-quality system in embodiment；
Fig. 3 is the response curve of the Three Degree Of Freedom time-varying system in embodiment.Wherein, Fig. 3 (A), Fig. 3 (B), figure
3 (C) represent mass m in Fig. 2 respectively1、m2、m3Response curve；
Fig. 4 is in embodiment under given arbitrary excitation, for Bayesian Estimation method and least-squares estimation
Method, using BIC criterion, AR exponent number and the exponent number of basic function.Wherein, Fig. 4 (A), Fig. 4 (B) are estimated using Bayes respectively
The exponent number of meter method and the least square estimation method AR exponent number and basic function；
Fig. 5 is for time-varying model, under 100 arbitrary excitations, using shellfish disclosed by the invention in embodiment
The system mode frequency that this method of estimation of leaf and the least square estimation method obtain, wherein, Fig. 5 (A), Fig. 5 (B) they are naTake 2, pa
15 are taken using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, is obtained
System mode frequency, Fig. 5 (C), Fig. 5 (D) are naTake 2, pa20 are taken to use Bayesian Estimation method disclosed by the invention and use
The least square estimation method, data length 8000, obtained system mode frequency, Fig. 5 (E), Fig. 5 (F) are naTake 2, paTake 25
Using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, what is obtained is
System modal frequency, Fig. 5 (G), Fig. 5 (H) are naTake 2, pa35 are taken using Bayesian Estimation method disclosed by the invention and using minimum
Two multiply method of estimation, and data length 8000, obtained system mode frequency, Fig. 5 (H), Fig. 5 (I) are naTake 2, paTake 60 uses
Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, obtained system mould
State frequency, Fig. 5 (J), Fig. 5 (K) are naTake 2, pa80 are taken using Bayesian Estimation method disclosed by the invention and using least square
Method of estimation, data length 8000, obtained system mode frequency；
Fig. 6 is for time-varying model, under 100 arbitrary excitations, using shellfish disclosed by the invention in embodiment
The system mode frequency that this method of estimation of leaf and the least square estimation method obtain, wherein, Fig. 6 (A), Fig. 6 (B) they are naTake 6, pa
15 are taken using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, is obtained
System mode frequency, Fig. 6 (C), Fig. 6 (D) are naTake 6, pa25 are taken to use Bayesian Estimation method disclosed by the invention and use
The least square estimation method, data length 8000, obtained system mode frequency, Fig. 6 (E), Fig. 6 (F) are naTake 6, paTake 35
Using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, what is obtained is
System modal frequency；
Fig. 7 is for time-varying model, under 100 arbitrary excitations, using shellfish disclosed by the invention in embodiment
The system mode frequency that this method of estimation of leaf and the least square estimation method obtain, wherein, Fig. 7 (A), Fig. 7 (B) they are naTake 3, pa
25 are taken using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, is obtained
System mode frequency, Fig. 7 (C), Fig. 7 (D) are naTake 4, pa25 are taken to use Bayesian Estimation method disclosed by the invention and use
The least square estimation method, data length 8000, obtained system mode frequency, Fig. 7 (E), Fig. 7 (F) are naTake 5, paTake 25
Using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 8000, what is obtained is
System modal frequency；
Fig. 8 is for time-varying model, under 100 arbitrary excitations, using shellfish disclosed by the invention in embodiment
The system mode frequency that this method of estimation of leaf and the least square estimation method obtain, wherein, Fig. 8 (A), Fig. 8 (B) they are naTake 2, pa
25 are taken using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 1000, is obtained
System mode frequency, Fig. 8 (C), Fig. 8 (F) are naTake 4, pa25 are taken to use Bayesian Estimation method disclosed by the invention and use
The least square estimation method, data length 1000, obtained system mode frequency, Fig. 8 (E), Fig. 8 (F) are naTake 6, paTake 25
Using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 1000, what is obtained is
System modal frequency；
Fig. 9 is for time-varying model, under 100 arbitrary excitations, using shellfish disclosed by the invention in embodiment
The system mode frequency that this method of estimation of leaf and the least square estimation method obtain, wherein, Fig. 9 (A), Fig. 9 (B) they are naTake 2, pa
30 are taken using Bayesian Estimation method disclosed by the invention and using the least square estimation method, data length 1000, is obtained
System mode frequency, Fig. 9 (C), Fig. 9 (D) are naTake 2, pa25 are taken to use Bayesian Estimation method disclosed by the invention and use
The least square estimation method, data length 1000, obtained system mode frequency；
Embodiment
In order to which objects and advantages of the present invention are better described, during below by Three Degree Of Freedom under an arbitrary excitation
Structure changes carry out dynamic analysis, and the present invention is made and explained in detail.
Embodiment 1：
Three Degree Of Freedom spring-dampers-the quality system of the present embodiment, as shown in Figure 2.The Three Degree Of Freedom bullet of the present embodiment
Spring-damper-quality system includes three mass m1、m2、m3, four damper c1、c2、c3、c4, four spring k1(t)、k2
(t)、k3(t)、k4(t), wherein, three masses and four dampers are permanent, are not changed over time, and four springs
Rigidity change over time.
Rigidity k in Three Degree Of Freedom time-varying systemi(t) change over time shown in relation such as formula (26)：
Each variable-value of system with 3 degrees of freedom is as shown in table 1.
The Three Degree Of Freedom time-varying system parameter of table 1
The dynamic control equation of system with 3 degrees of freedom shown in Fig. 2 is as shown in formula：
In formula, M, C are respectively mass matrix and damping matrix, and K (t) is the stiffness matrix of time-varying, and x is the response of system, f
(t) it is dynamic excitation.M, C, K (t) are as shown in formula (28)：
Gauss white-noise excitations are applied to three masses respectively, obtain the response of three masses in system.System is rung
Variable step Fourth order Runge-Kutta should be used to solve.Three degree of freedom displacement is used to recognize response signal samples used, will
For the solving result of runge kutta method with f=16Hz resamplings, the record time is 500s (t ∈ [0,500]), signal length N=
8000, output channel N0=3, the dynamic respond curve of three degree of freedom is as shown in Figure 3.
As shown in figure 1, a kind of only output Time variable structure modal parameter Bayesian Estimation method disclosed in the present embodiment, including
Following steps：
Step 1：Stochastic linear Time variable structure is represented using functional sequence multi -components time correlation autoregression model FS-VTAR
System.
To this embodiment, functional sequence multi -components time correlation autoregression model FS-VTAR expression shape is drawn
Shown in formula such as formula (29)：
Y=Φ θ+ε (29)
In formula,
Step 2：Using the model parameter θ and covariance matrix ∑ of FS-VTAR models in bayes method estimating step 1
Posterior distrbutionp.
Calculate the response Y at all moment, i.e. joint probability density function：
Shown in the condition prior distribution such as formula (34) for calculating θ：
In formula, m0It is parameter matrix θ expected matrix, M0It is covariance matrix,Expression average is a, association side
Poor matrix is the 3 of b2napaTie up normal distribution,Represent Kronecker product.For positive definite matrix ∑, prior distribution is taken as inverse
Wishart distribution, as shown in formula (35)：
∑~IW3×3(b0,Ψ0) (35)
In formula, IW3×3(b0,Ψ0) represent that it is b to obey the free degree0, Scale Matrixes Ψ0The N of distributionoRank positive definite real matrix.
In the case where θ and ∑ are given, all time of day response Y probability density function, i.e. likelihood function are calculated, such as
Shown in formula (36)：
In formula, the mark of tr () representing matrix, exp () represents natural Exponents, and p (x | y) represents the bar of x when known to y
Part probability density function, | | represent determinant.
θ Posterior distrbutionp is calculated：
With the Posterior distrbutionp of ∑：
In formula, parameter b and Ψ are respectively：
Step 3：The model order of functional sequence multi -components time correlation autoregression model is determined according to bayesian information criterion
Number and functional space exponent number, i.e. naAnd pa
According to bayesian information criterion BIC, as model order naWith functional space exponent number paWhen being optimal, BIC value takes
Minimum value.In practical operation, the response signal once obtained is taken, builds the complete functional sequence of vectors as shown in formula (39.)
Time-varying autoregressive, and take different naValue and paValue is modeled, to each group of naValue and paValue calculates BIC value, and will
It is plotted on same figure, is chosen BIC and is reached corresponding n during minimum valueaValue and paValue, it is selected exponent number.Fig. 4 is illustrated
Under given arbitrary excitation, for Bayesian Estimation and least-squares estimation, responded using different caused by BIC criterion.Functional
When sequence multi -components time correlation autoregression model FS-VTAR exponent numbers are optimal, BIC value takes minimum.Therefore, BIC values take
Corresponding functional sequence of vectors time-varying autoregressive FS-VTAR exponent numbers are n when minimumaValue.As shown in Figure 4, na=2 be pair
In the optimal selection of AR models.Corresponding to na=2, paBetween desirable 20 to 25.
Step 4：Time variable structure modal parameter and its credibility interval are asked using time freezing method.
Step 4.1：The calculating time freezes Time variable structure modal parameter.
System pole is calculated：
(D[k]-λr[k]I)Vr[k]=0 (40)
In formula, λr[k] andBe respectively k moment D [k] r-th of characteristic value and feature to
Amount, D [k] can be obtained by AR coefficient matrixes：
" time freezes " the modal frequency f of system is calculatedr[t] and dampingratioζr[t]：
Realize the Modal Parameter Identification for only output linearity Time variable structure.
Step 4.2：The credibility interval of Time variable structure modal parameter is obtained using Monte Carlo method
FromSampling, carry out Monte Carlo simulation.Divide for convenience of from probability
ClothSampling, construct random matrixOr
In formula, matrix Z each element is satisfied by unit normal distribution, CMAnd C∑It is Cholesky points of M and ∑ respectively
Solution.
∑ is obtained so as to build W by the expectation E (∑) of covariance matrix ∑：
Can be fromThe sample W that there is same distribution with parameter matrix θ is obtained, even θ=W.
Obtained parameter matrix θ is sampled every time, can obtain time frozen structure modal parameter, and this process is designated as a Monte Carlo mould
Intend.This step is repeated, the multiple Monte Carlo simulation of Time variable structure modal parameter Bayesian Estimation method is carried out, afterwards using point
The method of estimation can obtain estimating the credibility interval of modal parameter, obtain estimating not in the identification of linear time-varying modal parameters
Know the uncertainty of parameter, and then have and have the advantage that：The robustness for crossing estimation can be improved；Reduce and model structure is become
The susceptibility of change, the shorter response data of more accurate processing.The average value of Time variable structure modal parameter is repeatedly to cover spy
The average value of Carlow analog result, variance are the variance of multiple Monte Carlo simulation result.
Under 100 arbitrary excitations, using Bayesian Estimation method disclosed by the invention and the least square estimation method, na
2,3,4,5,6, p are taken respectivelya15,20,25,35,60,80 are taken respectively, calculate time-varying system modal parameter.As a result such as figure (5) (6)
(7) shown in.From scheming (5), when AR models take optimal, i.e. naWhen=2, basic function exponent number takes 15,20,25,35,60,80, most
Small square law and bayes method disclosed by the invention can realize good parameter Estimation, when large, disclosed by the invention
Bayes method estimated result is accurate, and least square method shows slight crossing and estimated.From scheming (6) (7), when AR rank
Number is 6, and basic function exponent number takes 15, when 25,35, the estimated result and n of bayes method disclosed by the inventionaVery phase when=2
Seemingly, only some magnitude frequencies are located at specific limit frequency band.And the identification result of least square method is with naIncrease become
Worse and worse, magnitude frequency is distributed among whole frequency band.From scheming (5) (6) (7), bayes method disclosed by the invention
It can be good at controlling estimation phenomenon to occur compared to existing least square method, have well in Structural Dynamics field
Application prospect.
Data length is changed into 1000 from 8000, exiting form is 100 arbitrary excitations, naAnd paTake different value.As a result exist
Scheme to present in (8) and (9).Work as naAnd paWhen taking optimal value, the estimated result obtained using least square method is not so good as using the present invention
Disclosed Bayesian Estimation method is preferable, but can still receive.And work as naAnd paWhen becoming big, obtained using least square method
Estimated result is very poor.And use the result and obtained when data volume is 8000 that Bayesian Estimation method disclosed by the invention obtains
Result it is almost consistent, the results showed that Bayesian Estimation method disclosed by the invention can be more compared to existing least square method
Change for the shorter response data of accurately processing and for model structure shows good robustness, in structural dynamic
Field has a good application prospect.
A kind of structure for only exporting the identification of Time variable structure modal parameter Bayesian Estimation method according to disclosed in the present embodiment
Modal parameter, the vibration frequency range of engineering structure can be obtained, detect whether to meet the requirement such as engineering structure standard and vibration isolation, energy
Enough instruct the design of engineering structure.In addition, obtained modal parameters can also be health monitoring, the structure failure of Time variable structure
Diagnosis, the application of structural vibration control etc. provide strong support, are with a wide range of applications and benefit.
Above-described specific descriptions, the purpose, technical scheme and beneficial effect of invention are carried out further specifically
It is bright, the specific embodiment that the foregoing is only the present invention is should be understood that, for explaining the present invention, is not used to limit this
The protection domain of invention, within the spirit and principles of the invention, any modification, equivalent substitution and improvements done etc. all should
Within protection scope of the present invention.
Claims (6)
1. a kind of only export Time variable structure modal parameter Bayesian Estimation method, it is characterised in that：Comprise the following steps,
Step 1：Linear time-varying structural system is represented using functional sequence multi -components time correlation autoregression model FS-VTAR；
Step 1 concrete methods of realizing, comprises the following steps：
Step 1.1：Represented with the form of time discrete shown in multi -components time correlation autoregression model VTAR such as formulas (1)：
<mrow>
<mi>x</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<mo>-</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
</munderover>
<msub>
<mi>A</mi>
<mi>i</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mi>x</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>-</mo>
<mi>i</mi>
<mo>&rsqb;</mo>
<mo>+</mo>
<mi>e</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>,</mo>
<mi>e</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>-</mo>
<mi>N</mi>
<mi>I</mi>
<mi>D</mi>
<mrow>
<mo>(</mo>
<mn>0</mn>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>(</mo>
<mi>k</mi>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula,It is NoThe response vector signal of output channel, e [k] is can not as caused by measurement and model approximation
The non-stationary residual error of observation, e [k] average is zero, and time-varying covariance matrix isnaFor VTAR model orders,
Ai[k] is the i-th rank time-varying AR coefficient matrixes, and k represents k-th of discrete time point；
Step 1.2：By time-varying AR factor projections into a series of function space represented by given basic function sequences；
Shown in given function space such as formula (2)：
In formula, p be function space dimension, fjFor jth rank basis function vector；
By Ai[k] deploys, that is, obtains function spaceThe functional sequence of vectors time-varying autoregressive (FS-VTAR) of expression, such as
Shown in formula (3)：
<mrow>
<mi>x</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<mo>-</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
</munderover>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>p</mi>
<mi>a</mi>
</msub>
</munderover>
<msub>
<mi>A</mi>
<mrow>
<mi>i</mi>
<mi>j</mi>
</mrow>
</msub>
<msub>
<mi>f</mi>
<mi>j</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mi>x</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>-</mo>
<mi>i</mi>
<mo>&rsqb;</mo>
<mo>+</mo>
<mi>e</mi>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
Shown in x [k] component form such as formula (4)：
<mrow>
<msub>
<mi>x</mi>
<mi>l</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<mo>-</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
</munderover>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>p</mi>
<mi>a</mi>
</msub>
</munderover>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>m</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</munderover>
<msub>
<mi>a</mi>
<mrow>
<mi>i</mi>
<mi>j</mi>
<mo>,</mo>
<mi>l</mi>
<mi>m</mi>
</mrow>
</msub>
<msub>
<mi>f</mi>
<mi>j</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<msub>
<mi>x</mi>
<mi>m</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>-</mo>
<mi>i</mi>
<mo>&rsqb;</mo>
<mo>+</mo>
<msub>
<mi>e</mi>
<mi>l</mi>
</msub>
<mo>&lsqb;</mo>
<mi>k</mi>
<mo>&rsqb;</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula, aij,lmFor matrix AijElement, l represent l-th of output channel；
Step 1.3：Shown in the expression-form such as formula (5) for drawing functional sequence multi -components time correlation autoregression model FS-VTAR：
Y=Φ θ+ε (5)
In formula,
Wherein, symbol " T " represents transposition computing；
So far, complete to represent linear time-varying structural system using functional sequence multi -components time correlation autoregression model FS-VTAR；
Step 2：Using the model parameter θ of FS-VTAR models in bayes method estimating step 1 and the posteriority of covariance matrix ∑
Distribution；
Step 3：According to bayesian information criterion determine functional sequence multi -components time correlation autoregression model model order and
N in functional space exponent number, i.e. formula (4)aAnd pa；
Step 4：Time variable structure modal parameter and credibility interval are asked using time freezing method.
A kind of 2. only output Time variable structure modal parameter Bayesian Estimation method as claimed in claim 1, it is characterised in that：Also
Including step 5,
The modal parameter that applying step 1 obtains to step 4 identification instructs the Structural analysis and design in Structural Dynamics field.
3. a kind of only output Time variable structure modal parameter Bayesian Estimation method as claimed in claim 1 or 2, its feature exist
In：Step 5 concrete methods of realizing is,
The modal parameter obtained according to step 1 to step 4 identification, can obtain the vibration frequency range of engineering structure, detect whether
Meet the requirement such as engineering structure standard and vibration isolation, the design of engineering structure can be instructed；In addition, obtained modal parameters are also
Strong support can be provided for the application of the health monitoring, structure failure diagnosis, structural vibration control of Time variable structure etc., had
Have wide practical use and benefit.
A kind of 4. only output Time variable structure modal parameter Bayesian Estimation method as claimed in claim 3, it is characterised in that：Step
Rapid 2 concrete methods of realizing comprises the following steps,
Step 2.1：Obtain the response Y at all moment, i.e. formula (5) left end item joint probability density function；
Covariance matrix ∑ (k) in formula (1) is definite value, i.e. ∑ (t)=∑；In the viewpoint of Frequency school, unknown parameter θ and
∑ is definite value, and bayes method thinks that unknown parameter θ and ∑ are stochastic variable；Therefore, the distribution for responding Y is with parameter
θ and the probability distribution that ∑ is condition, its probability density function is represented by a multivariate normal distributions, as shown in formula (9)：
<mrow>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mrow>
<mo>(</mo>
<mi>Y</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>~</mo>
<msub>
<mi>N</mi>
<mrow>
<mo>(</mo>
<mi>N</mi>
<mo>-</mo>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
<mo>)</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mo>(</mo>
<mrow>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>9</mn>
<mo>)</mo>
</mrow>
</mrow>
Step 2.2：Obtain the prior distribution of parameter matrix θ and covariance matrix ∑；
Parameter matrix θ and covariance matrix ∑ is unknowable, based on linear model standard conjugate gradient descent method theory, in formula
(1) when the e [k] in is given as normal distribution, under ∑ known case, shown in θ condition prior distribution such as formula (10)：
<mrow>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>&Sigma;</mi>
<mo>~</mo>
<msub>
<mi>N</mi>
<mrow>
<msubsup>
<mi>N</mi>
<mi>o</mi>
<mn>2</mn>
</msubsup>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
<msub>
<mi>p</mi>
<mi>a</mi>
</msub>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mo>(</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>&CircleTimes;</mo>
<msub>
<mi>M</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>10</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula,It is parameter matrix θ expected matrix,It is covariance matrix,Table
Show average be a, covariance matrix be b'sTie up normal distribution,Represent Kronecker product；For positive definite matrix ∑,
Prior distribution is taken as inverse Wishart distribution, as shown in formula (11)：
<mrow>
<mi>&Sigma;</mi>
<mo>~</mo>
<msub>
<mi>IW</mi>
<mrow>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>&times;</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>,</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>11</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula,Represent that it is b to obey the free degree0, Scale Matrixes Ψ0The N of distributionoRank positive definite real matrix；
Step 2.3：Obtain the likelihood function and priori probability density distribution function in bayes method；
From formula (9), in the case where θ and ∑ are given, all time of day response Y probability density function, i.e. likelihood letter
Number, as shown in formula (12)：
<mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>r</mi>
<mo>|</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>NN</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mi>&Sigma;</mi>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mi>N</mi>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mi>&Sigma;</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>12</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula, the mark of tr () representing matrix, exp () expression natural Exponents, p (x | y) represent that x condition is general when known to y
Rate density function, | | represent determinant；From formula (10) and multivariate normal distributions equation, the θ in the case where ∑ is given
Priori probability density function p (θ | ∑) it is expressed as：
<mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>|</mo>
<mi>&Sigma;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mi>&Sigma;</mi>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mi>&Sigma;</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<mi>m</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>13</mn>
<mo>)</mo>
</mrow>
</mrow>
Based on formula (11), the priori probability density function p (∑) of ∑ is expressed as：
<mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&Sigma;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mrow>
<mo>|</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<msup>
<mo>|</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
</msup>
</mrow>
<mrow>
<msup>
<mn>2</mn>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<msub>
<mi>&Gamma;</mi>
<mi>p</mi>
</msub>
<mrow>
<mo>(</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>|</mo>
<mi>&Sigma;</mi>
<msup>
<mo>|</mo>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mi>&Sigma;</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>14</mn>
<mo>)</mo>
</mrow>
</mrow>
Then from formula (13) and formula (14), the joint priori probability density function p (θ, ∑) of θ and ∑ is：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mo>&Sigma;</mo>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>|</mo>
<mo>&Sigma;</mo>
<mo>)</mo>
</mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mo>&Sigma;</mo>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>&times;</mo>
<mfrac>
<mrow>
<mo>|</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<msup>
<mo>|</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
</msup>
</mrow>
<mrow>
<msup>
<mn>2</mn>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<msub>
<mi>&Gamma;</mi>
<mi>p</mi>
</msub>
<mrow>
<mo>(</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>15</mn>
<mo>)</mo>
</mrow>
</mrow>
Step 2.4：The parameter matrix θ and the Posterior distrbutionp of covariance matrix ∑ that calculating bayes method obtains；
Because response Y is unrelated with θ and ∑, p (Y | θ, ∑) marginal likelihood is constant, in the case of Y determinations, after θ and ∑
Probability density function is tested to be expressed as：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mo>&Sigma;</mo>
<mo>|</mo>
<mi>Y</mi>
<mo>)</mo>
</mrow>
<mo>&Proportional;</mo>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>Y</mi>
<mo>|</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mo>&Sigma;</mo>
<mo>)</mo>
</mrow>
<mi>p</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>,</mo>
<mo>&Sigma;</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>NN</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mi>N</mi>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>&times;</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>&times;</mo>
<mfrac>
<mrow>
<mo>|</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<msup>
<mo>|</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
</msup>
</mrow>
<mrow>
<msup>
<mn>2</mn>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<msub>
<mi>&Gamma;</mi>
<mi>p</mi>
</msub>
<mrow>
<mo>(</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>NN</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mfrac>
<mrow>
<mo>|</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<msup>
<mo>|</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
</msup>
</mrow>
<mrow>
<msup>
<mn>2</mn>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<msub>
<mi>&Gamma;</mi>
<mi>p</mi>
</msub>
<mrow>
<mo>(</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<mi>N</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>&times;</mo>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
</mrow>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>Y</mi>
<mo>-</mo>
<mi>&Phi;</mi>
<mi>&theta;</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<mfrac>
<mrow>
<mo>|</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<msup>
<mo>|</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
</msup>
</mrow>
<mrow>
<msup>
<mn>2</mn>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</msup>
<msub>
<mi>&Gamma;</mi>
<mi>p</mi>
</msub>
<mrow>
<mo>(</mo>
<mfrac>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mn>2</mn>
</mfrac>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mrow>
<mo>(</mo>
<mrow>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msubsup>
<mi>m</mi>
<mn>0</mn>
<mi>T</mi>
</msubsup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msup>
<mi>Y</mi>
<mi>T</mi>
</msup>
<mi>Y</mi>
<mo>-</mo>
<msup>
<mi>m</mi>
<mi>T</mi>
</msup>
<msup>
<mi>M</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mi>m</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>&times;</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mi>&pi;</mi>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>NN</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mo>|</mo>
<mo>&Sigma;</mo>
<msup>
<mo>|</mo>
<mrow>
<mo>-</mo>
<mfrac>
<mrow>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>+</mo>
<mi>N</mi>
<mo>+</mo>
<mn>2</mn>
</mrow>
<mn>2</mn>
</mfrac>
</mrow>
</msup>
<mi>exp</mi>
<mrow>
<mo>(</mo>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mi>t</mi>
<mi>r</mi>
<mo>(</mo>
<mrow>
<msup>
<mo>&Sigma;</mo>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<mi>m</mi>
</mrow>
<mo>)</mo>
</mrow>
<mi>T</mi>
</msup>
<msup>
<mi>M</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mrow>
<mo>(</mo>
<mrow>
<mi>&theta;</mi>
<mo>-</mo>
<mi>m</mi>
</mrow>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>16</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula,
<mrow>
<msup>
<mi>M</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mo>=</mo>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>+</mo>
<msup>
<mi>&Phi;</mi>
<mi>T</mi>
</msup>
<mi>&Phi;</mi>
<mo>,</mo>
<mi>m</mi>
<mo>=</mo>
<mi>M</mi>
<mrow>
<mo>(</mo>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msup>
<mi>&Phi;</mi>
<mi>T</mi>
</msup>
<mi>Y</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>17</mn>
<mo>)</mo>
</mrow>
</mrow>
Contrast (16) and formula (15), obtains θ Posterior distrbutionp：
<mrow>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mrow>
<mo>(</mo>
<mi>&theta;</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mi>r</mi>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>~</mo>
<msub>
<mi>N</mi>
<mrow>
<msubsup>
<mi>N</mi>
<mi>o</mi>
<mn>2</mn>
</msubsup>
<msub>
<mi>n</mi>
<mi>a</mi>
</msub>
<msub>
<mi>p</mi>
<mi>a</mi>
</msub>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>v</mi>
<mi>e</mi>
<mi>c</mi>
<mo>(</mo>
<mi>m</mi>
<mo>)</mo>
<mo>,</mo>
<mi>&Sigma;</mi>
<mo>&CircleTimes;</mo>
<mi>M</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>18</mn>
<mo>)</mo>
</mrow>
</mrow>
With the Posterior distrbutionp of ∑：
<mrow>
<mi>&Sigma;</mi>
<mo>|</mo>
<mi>Y</mi>
<mo>,</mo>
<mo>~</mo>
<msub>
<mi>IW</mi>
<mrow>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>&times;</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>b</mi>
<mo>,</mo>
<mi>&Psi;</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>19</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula, parameter b and Ψ are respectively
<mrow>
<mi>&Psi;</mi>
<mo>=</mo>
<msub>
<mi>&Psi;</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msubsup>
<mi>m</mi>
<mn>0</mn>
<mi>T</mi>
</msubsup>
<msubsup>
<mi>M</mi>
<mn>0</mn>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<msub>
<mi>m</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<msup>
<mi>Y</mi>
<mi>T</mi>
</msup>
<mi>Y</mi>
<mo>-</mo>
<msup>
<mi>m</mi>
<mi>T</mi>
</msup>
<mi>M</mi>
<mi>m</mi>
<mo>,</mo>
<mi>b</mi>
<mo>=</mo>
<msub>
<mi>b</mi>
<mn>0</mn>
</msub>
<mo>+</mo>
<mi>N</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>20</mn>
<mo>)</mo>
</mrow>
</mrow>
A kind of 5. only output Time variable structure modal parameter Bayesian Estimation method as claimed in claim 4, it is characterised in that：Step
Rapid 3 concrete methods of realizing is as follows,
According to bayesian information criterion BIC, as model order naWith functional space exponent number paWhen being optimal, BIC value takes minimum
Value；In practical operation, take the response signal once obtained, build complete functional sequence of vectors time-varying as shown in formula (6) from
Regression model, and take different naValue and paValue is modeled, to each group of naValue and paValue calculates BIC value, and is drawn
On same figure, choose BIC and reach corresponding n during minimum valueaValue and paValue, it is selected exponent number.
A kind of 6. only output Time variable structure modal parameter Bayesian Estimation method as claimed in claim 5, it is characterised in that：Step
Rapid 4 concrete methods of realizing comprises the following steps,
Step 4.1：The calculating time freezes Time variable structure modal parameter；
Moment point k is freezed using the method for " time freezes ", by VTAR coefficients Ai[k] obtains modal parameter；Using adjoint matrix by
Generalized eigenvalue problem obtains system pole：
(D[k]-λr[k]I)Vr[k]=0 (21)
In formula, λr[k] andIt is r-th of characteristic value and characteristic vector at k moment D [k] respectively, by
AR coefficient matrixes obtain D [k]：
" time freezes " the modal frequency f of system is calculatedr[t] and dampingratioζr[t]：
<mrow>
<msub>
<mi>f</mi>
<mi>r</mi>
</msub>
<mo>&lsqb;</mo>
<mi>t</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfrac>
<mrow>
<mn>2</mn>
<mi>&pi;</mi>
<mo>|</mo>
<msub>
<mi>ln&lambda;</mi>
<mi>r</mi>
</msub>
<mo>&lsqb;</mo>
<mi>t</mi>
<mo>&rsqb;</mo>
<mo>|</mo>
</mrow>
<msub>
<mi>T</mi>
<mi>s</mi>
</msub>
</mfrac>
<mo>,</mo>
<msub>
<mi>&zeta;</mi>
<mi>r</mi>
</msub>
<mo>&lsqb;</mo>
<mi>t</mi>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfrac>
<mrow>
<mo>-</mo>
<mi>Re</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>ln&lambda;</mi>
<mi>r</mi>
</msub>
<mo>&lsqb;</mo>
<mi>t</mi>
<mo>&rsqb;</mo>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mo>|</mo>
<msub>
<mi>ln&lambda;</mi>
<mi>r</mi>
</msub>
<mo>&lsqb;</mo>
<mi>t</mi>
<mo>&rsqb;</mo>
<mo>|</mo>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>23</mn>
<mo>)</mo>
</mrow>
</mrow>
Realize the Modal Parameter Identification for only output linearity Time variable structure；
Step 4.2：The credibility interval of Time variable structure modal parameter is obtained using Monte Carlo method, obtains linear time-varying structural modal
The uncertainty of unknown parameter is estimated in parameter identification, tool has the advantage that：The robustness for crossing estimation can be improved；Reduce
To the susceptibility of model structure change, the shorter response data of more accurate processing；
FromSampling, carry out Monte Carlo simulation；For convenience from probability distributionSampling, random matrix is constructed by formula (24)Or
<mrow>
<mi>W</mi>
<mo>=</mo>
<mi>m</mi>
<mo>+</mo>
<msubsup>
<mi>C</mi>
<mi>M</mi>
<mrow>
<mo>-</mo>
<mi>T</mi>
</mrow>
</msubsup>
<msubsup>
<mi>ZC</mi>
<mi>&Sigma;</mi>
<mi>T</mi>
</msubsup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>24</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula, matrix Z each element is satisfied by unit normal distribution, CMAnd C∑It is the Cholesky decomposition of M and ∑ respectively；
∑ is obtained so as to build W by the expectation E (∑) of covariance matrix ∑, according to inverse Wishart distribution, i.e. formula (25)：
<mrow>
<mi>E</mi>
<mrow>
<mo>(</mo>
<mi>&Sigma;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mi>&Psi;</mi>
<mrow>
<mi>b</mi>
<mo>-</mo>
<msub>
<mi>N</mi>
<mi>o</mi>
</msub>
<mo>-</mo>
<mn>1</mn>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>25</mn>
<mo>)</mo>
</mrow>
</mrow>
By formula (24) and (25), fromThe sample W that there is same distribution with parameter matrix θ is obtained, i.e.,
Make θ=W；Obtained parameter matrix θ is sampled every time, time frozen structure modal parameter, this process can be obtained according to formula (23)
It is designated as single Monte Carlo simulation；Single Monte Carlo simulation step is repeated, carries out Time variable structure modal parameter Bayesian Estimation
The multiple Monte Carlo simulation of method, obtain estimating the credibility interval of modal parameter using the method for point estimation afterwards, obtain line
Property Time variable structure Modal Parameter Identification in estimate unknown parameter uncertainty, and then have have the advantage that：Can improve for
Cross the robustness of estimation；Reduce the susceptibility to model structure change, the shorter response data of more accurate processing；Time-varying knot
The average value of structure modal parameter is the average value of multiple Monte Carlo simulation result, and variance is multiple Monte Carlo simulation knot
The variance of fruit.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710875350.6A CN107844627A (en) | 2017-09-25 | 2017-09-25 | It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710875350.6A CN107844627A (en) | 2017-09-25 | 2017-09-25 | It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN107844627A true CN107844627A (en) | 2018-03-27 |
Family
ID=61661416
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710875350.6A Pending CN107844627A (en) | 2017-09-25 | 2017-09-25 | It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107844627A (en) |
Cited By (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109598027A (en) * | 2018-11-08 | 2019-04-09 | 合肥工业大学 | A kind of algorithm based on frequency response function correcting principle model parameter |
| CN109635389A (en) * | 2018-11-29 | 2019-04-16 | 中国航空工业集团公司沈阳飞机设计研究所 | A kind of electric steering engine stiffness test data processing method |
| CN110489891A (en) * | 2019-08-23 | 2019-11-22 | 江南大学 | A kind of industrial process time-varying uncertainty method based on more born of the same parents' space filterings |
| CN111273297A (en) * | 2019-03-06 | 2020-06-12 | 哈尔滨工程大学 | Horizontal array moving target depth estimation method based on AR wave number spectrum |
| CN111859250A (en) * | 2020-06-19 | 2020-10-30 | 哈尔滨工业大学(深圳)(哈尔滨工业大学深圳科技创新研究院) | Method, system and storage medium for evaluating accuracy of rapid modal parameters under environment and forced excitation |
| CN111884859A (en) * | 2020-07-30 | 2020-11-03 | 国网重庆市电力公司电力科学研究院 | Network fault diagnosis method and device and readable storage medium |
| CN112154464A (en) * | 2018-06-19 | 2020-12-29 | 株式会社岛津制作所 | Parameter search method, parameter search device, and program for parameter search |
| CN112837816A (en) * | 2021-02-09 | 2021-05-25 | 清华大学 | Physiological state prediction method, computer device, and storage medium |
| CN113822354A (en) * | 2021-09-17 | 2021-12-21 | 合肥工业大学 | Micro-nano probe dynamic characteristic compensation method based on Bayesian inverse calculus modeling |
| CN114048548A (en) * | 2021-11-18 | 2022-02-15 | 汉思科特(盐城)减震技术有限公司 | Method for estimating dynamic load of air spring of heavy-duty vehicle |
| WO2022105104A1 (en) * | 2020-11-18 | 2022-05-27 | 南通大学 | Multi-innovation recursive bayesian algorithm-based battery model parameter identification method |
-
2017
- 2017-09-25 CN CN201710875350.6A patent/CN107844627A/en active Pending
Cited By (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112154464B (en) * | 2018-06-19 | 2024-01-02 | 株式会社岛津制作所 | Parameter searching method, parameter searching device, and parameter searching program |
| CN112154464A (en) * | 2018-06-19 | 2020-12-29 | 株式会社岛津制作所 | Parameter search method, parameter search device, and program for parameter search |
| CN109598027A (en) * | 2018-11-08 | 2019-04-09 | 合肥工业大学 | A kind of algorithm based on frequency response function correcting principle model parameter |
| CN109635389A (en) * | 2018-11-29 | 2019-04-16 | 中国航空工业集团公司沈阳飞机设计研究所 | A kind of electric steering engine stiffness test data processing method |
| CN109635389B (en) * | 2018-11-29 | 2022-12-20 | 中国航空工业集团公司沈阳飞机设计研究所 | Rigidity test data processing method for electric steering engine |
| CN111273297A (en) * | 2019-03-06 | 2020-06-12 | 哈尔滨工程大学 | Horizontal array moving target depth estimation method based on AR wave number spectrum |
| CN110489891B (en) * | 2019-08-23 | 2020-11-17 | 江南大学 | Industrial process time-varying parameter estimation method based on multi-cell spatial filtering |
| CN110489891A (en) * | 2019-08-23 | 2019-11-22 | 江南大学 | A kind of industrial process time-varying uncertainty method based on more born of the same parents' space filterings |
| CN111859250A (en) * | 2020-06-19 | 2020-10-30 | 哈尔滨工业大学(深圳)(哈尔滨工业大学深圳科技创新研究院) | Method, system and storage medium for evaluating accuracy of rapid modal parameters under environment and forced excitation |
| CN111884859A (en) * | 2020-07-30 | 2020-11-03 | 国网重庆市电力公司电力科学研究院 | Network fault diagnosis method and device and readable storage medium |
| WO2022105104A1 (en) * | 2020-11-18 | 2022-05-27 | 南通大学 | Multi-innovation recursive bayesian algorithm-based battery model parameter identification method |
| CN112837816A (en) * | 2021-02-09 | 2021-05-25 | 清华大学 | Physiological state prediction method, computer device, and storage medium |
| CN113822354A (en) * | 2021-09-17 | 2021-12-21 | 合肥工业大学 | Micro-nano probe dynamic characteristic compensation method based on Bayesian inverse calculus modeling |
| CN113822354B (en) * | 2021-09-17 | 2022-12-06 | 合肥工业大学 | Micro-nano probe dynamic characteristic compensation method based on Bayesian inverse calculus modeling |
| CN114048548A (en) * | 2021-11-18 | 2022-02-15 | 汉思科特(盐城)减震技术有限公司 | Method for estimating dynamic load of air spring of heavy-duty vehicle |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107844627A (en) | It is a kind of only to export Time variable structure modal parameter Bayesian Estimation method | |
| CN106354695B (en) | One kind only exporting linear Time variable structure Modal Parameters Identification | |
| CN107590317B (en) | Generator dynamic estimation method considering model parameter uncertainty | |
| CN103900610B (en) | MEMS gyro random error Forecasting Methodology based on Lycoperdon polymorphum Vitt wavelet neural network | |
| CN107133195A (en) | A kind of Methodology of The Determination of The Order of Model of engineering structure modal idenlification | |
| CN106897717A (en) | Bayesian model modification method under multiple test based on environmental excitation data | |
| CN101916241B (en) | Method for identifying time-varying structure modal frequency based on time frequency distribution map | |
| CN105142177A (en) | Complex neural network channel prediction method | |
| CN103246890B (en) | Modal Parameters Identification based on multi-input multi-output signal noise reduction | |
| CN110705041B (en) | EASI-based linear structure working modal parameter identification method | |
| CN112086100B (en) | Quantization error entropy based urban noise identification method of multilayer random neural network | |
| CN105043384A (en) | Modeling method of gyroscopic random noise ARMA model based on robust Kalman wave filtering | |
| CN104463347A (en) | Method for electronic product degenerate state trend prediction with singular signals | |
| CN107729845A (en) | A kind of frequency respond noise-reduction method decomposed based on sub-space feature value | |
| CN104795063A (en) | Acoustic model building method based on nonlinear manifold structure of acoustic space | |
| Cerone et al. | Set-membership errors-in-variables identification of MIMO linear systems | |
| CN111274630A (en) | Physical mode extraction method for engineering structure flexibility recognition | |
| Mai et al. | Surrogate modelling for stochastic dynamical systems by combining NARX models and polynomial chaos expansions | |
| CN109254321B (en) | Method for identifying rapid Bayesian modal parameters under seismic excitation | |
| CN113721458A (en) | Nonlinear system fault estimation observation method and device with disturbance | |
| Guillaume et al. | OMAX–a combined experimental-operational modal analysis approach | |
| CN113884914B (en) | Method for estimating charge state of evanescent parallel Kalman filtering power battery | |
| Shaikh et al. | Wavelet Decomposition Impacts on Traditional Forecasting Time Series Models. | |
| Seledkova et al. | Forecasting characteristics of time series to support managerial decision making process in production-And-economic systems | |
| CN111428342B (en) | Random dynamic load identification method based on frequency domain spectrum decomposition |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20180327 |
|
| WD01 | Invention patent application deemed withdrawn after publication |