Industrial object-oriented Gaussian approximation filtering method for nonlinear system

Abstract

The invention relates to the technical field of computers, and discloses an industrial object-oriented Gaussian approximation filtering method for a nonlinear system. Mathematical modeling is carried out on a physical system to be subjected to state estimation, and a nonlinear state space model corresponding to the physical system is obtained; matching the first Gaussian distribution of the system state at the time t to a given distribution through a moment matching operation based on the non-linear state-space model to determine a predictive estimate and a predictive estimate error of the system state at the time t; acquiring a system observation value at the moment t; iteratively updating the optimal estimation and the optimal estimation error of the system state at the moment t based on a natural gradient descent method and a Fisher information matrix associated with the second Gaussian distribution of the system state at the moment t, and obtaining the optimal estimation and the optimal estimation error at the moment t until an iteration stop condition is met; linearization errors can be avoided, and more accurate system state estimation of a nonlinear system is provided.

Description

Technical Field

[0001] The present disclosure relates to the field of computer technology, and in particular to a Gaussian approximate filtering method for a nonlinear system facing an industrial object. Background Art

[0002] State estimation is used to reconstruct the internal state of a system from the noisy observation data of the system in order to control, monitor or predict the system. State estimation of industrial objects has important application value in many fields of modern industry, such as automation control, process monitoring and fault detection. With the development of manufacturing industry towards intelligence and automation, accurate estimation of the state of industrial systems has become a key task to ensure the stability of system operation and optimize production processes. The complexity, nonlinearity and changeable environmental interference of industrial objects make the state estimation results of conventional filtering methods inaccurate. Therefore, it is an urgent problem to propose an efficient and stable filtering method for industrial objects.

[0003] Current state estimation methods are generally based on Bayesian filtering theory, which targets the current state distribution given the available measurement values ​​so far. The current state distribution, called the posterior distribution, is usually calculated through two iterative steps. The two-step iteration is: prediction step and update step. In the prediction step, the Chapman-Kolmogorov equation is used to predict the state distribution forward, and the transition model is used to describe the time evolution of the system state. Based on the predicted prior distribution, the update step applies Bayes' theorem and uses the observation model to update the prior, which describes the relationship between the noisy measurement and the true state. Specifically, the prediction step and the update step can be expressed by the following formula:

[0004] Prediction steps:

[0005] p(x t |y 1:t-1 )=∫p(x t |x t-1 )p(x t-1 |y 1:t-1 )dx t-1 ;

[0006] Update steps:

[0007]

[0008] Among them, p(x t |y 1:t-1 ) means that given all the observation data y from time step 1 to time step t-1 1:t-1 Under the condition that, at time step t, the system state x t The probability distribution of p(x t |xt-1 ) represents the transition model, that is, from the system state x t-1 Transfer to system state x t The probability of p(x t-1 |y 1:t-1 ) means that given all the observation data y from time step 1 to time step t-1 1:t Under the condition that, at time step t-1, the system state x t-1 The probability distribution of p(x t |y 1:t ) means that given all the observation data y from time step 1 to time step t 1:t Under the condition that, at time step t, the system state x t The probability distribution of p(y t |x t ) represents the observation model, that is, given the system state x t Under the condition that t The conditional probability distribution of .

[0009] However, industrial systems are often nonlinear, and the above Bayesian filter theory cannot be analytically calculated in the case where the transition probability and observation probability are nonlinear or non-Gaussian. In order to make Bayesian filtering computable in practical applications, it is usually assumed that the state distribution at each time step is Gaussian. Then, Kalman Filter (KF) is performed, that is, the Gaussian posterior probability is analytically calculated by recursively updating its mean and covariance.

[0010] When dealing with linear systems, Gaussian filters can show excellent performance; however, when dealing with nonlinear systems, traditional Gaussian filters (such as Extended Kalman Filter (EKF) or Unscented Kalman Filter (UKF)) rely on linearization technology. These technologies control the computational complexity by approximately linearizing the nonlinear system, but they will introduce large linearization errors, resulting in a decrease in the accuracy of state estimation, which seriously affects the reliability and operating efficiency of the entire system. Therefore, in nonlinear systems, how to accurately and quickly perform state estimation has become an important issue that needs to be solved urgently. Summary of the invention

[0011] In view of this, the present disclosure proposes a Gaussian approximate filtering method for nonlinear systems of industrial objects, which can solve the problem that the estimation results of traditional state estimation methods have large linearization errors.

[0012] According to one aspect of the present disclosure, a Gaussian approximation filtering method for a nonlinear system of an industrial object is provided, the method comprising:

[0013] Mathematically modeling the physical system for which state estimation is to be performed to obtain a nonlinear state space model corresponding to the physical system;

[0014] Based on the nonlinear state space model, a first Gaussian distribution of the system state at time t is matched with a given distribution through a moment matching operation to determine a predicted estimate and a predicted estimate error of the system state at time t; the given distribution takes the optimal estimate of the system state at time t-1 as a mean and the optimal estimate error of the system state at time t-1 as a covariance matrix; the optimal estimate and the optimal estimate error of the system state at time 0 are initialized as preset values; and t is a positive integer;

[0015] Get the system observation value at time t;

[0016] Based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, the optimal estimate and the optimal estimate error of the system state at time t are iteratively updated until the iteration stop condition is met, and the optimal estimate and the optimal estimate error at time t are obtained; wherein the second Gaussian distribution has the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix.

[0017] In a possible implementation, the iterative updating of the optimal estimate and the optimal estimate error of the system state at time t based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t includes:

[0018] Obtaining an optimal estimate of the system state at the time t and an iterative update model of the optimal estimate error; wherein the iterative update model is obtained by stacking the Gaussian parameters of the second Gaussian distribution into a single column vector, determining the derivative with respect to the column vector, and then combining the column vector, the derivative with respect to the column vector, the Fisher information matrix, and the expression of the natural gradient descent method;

[0019] The predicted estimate and predicted estimate error of the system state at time t, the optimal estimate and optimal estimate error of the i-th iteration at time t, and the system observation value at time t are input into the iterative update model to obtain the optimal estimate and optimal estimate error of the i+1-th iteration; wherein i is an integer increasing from 0, and the optimal estimate and optimal estimate error at the 0th iteration are set based on a preset method.

[0020] In a possible implementation, the iterative update model is expressed by the following formula:

[0021]

[0022] in, represents the optimal estimation error of the i+1th iteration The inverse of represents the optimal estimate of the i+1th iteration; represents the prediction estimation error of the system state at time t; represents the predicted estimate of the system state at time t; represents the optimal estimation error of the i-th iteration at time t The inverse of represents the optimal estimate of the i-th iteration at time t; represents the preset loss function; represents the system state x at time t t It follows a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; y t represents the system observation value at time t; Represents the expected value calculation function.

[0023] In a possible implementation, the preset method includes setting the optimal estimate and the optimal estimate error at the 0th iteration by Laplace approximation, and the preset method is expressed by the following formula:

[0024]

[0025] in, represents the optimal estimate at the 0th iteration at time t, represents the inverse of the optimal estimation error at the 0th iteration at time t; Represents the system state x estimated based on the maximum a posteriori probability at time t t The best estimate of represents the system state x at time t t It obeys the Gaussian distribution with the prediction estimate at time t as the mean and the prediction estimate error at time t as the covariance matrix; it represents the preset loss function; y t represents the system observation value at time t; exp represents the exponential operation function; Indicates the search for The system state x with the largest value t ; Indicates the system status Department, About system status x t The second derivative of .

[0026] In a possible implementation, the loss function is one of the following: an observation likelihood function, a Beta divergence-induced likelihood function, a pseudo-Huber loss function, and a weighted observation likelihood function.

[0027] In a possible implementation, the nonlinear state space model is used to perform a moment matching operation to match the first Gaussian distribution of the system state at time t with a given distribution to determine a predicted estimate and a predicted estimate error of the system state at time t, which is expressed by the following formula:

[0028]

[0029] in, represents the predicted estimate of the system state at time t; P t|t-1 represents the prediction estimation error of the system state at time t; represents the given distribution; f represents the transfer matrix in the nonlinear state space model; ξ t-1 A matrix representing the process noise at time t-1 in the nonlinear state-space model; represents the expected value calculation function; Represents the variance calculation function.

[0030] In a possible implementation, the nonlinear state space model is expressed by the following formula:

[0031] x t+1 =f(x t )+ξ t ;

[0032]

[0033] in, represents the system state of the physical system, represents n-dimensional Euclidean space; is the system observation value of the physical system, represents m-dimensional Euclidean space, n and m are positive integers; represents the transfer matrix; represents the observation matrix; ξ t represents the process noise of the physical system; represents the observation noise of the physical system, ξ t and They are all approximately Gaussian distributions with a mean of 0, that is:

[0034]

[0035] Among them, 0 n×1 Represents a matrix with n rows and 1 column of elements equal to 0, 0 m×1 represents a matrix with m rows and 1 column with elements equal to 0. is the covariance matrix of process noise, is the covariance matrix of the observation noise; ξ at different times tIndependent of each other, ξ at different times t Independent of each other, and ξ at any time t and are independent of each other and of the initial state x0 of the physical system; the initial state x0 follows a Gaussian distribution in, and are the preset initial values ​​of random variables and the initial values ​​of the uncertainty of the initial state respectively.

[0036] In a possible implementation, the iteration stop condition is expressed by the following formula:

[0037]

[0038] in, represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i+1th iteration as the mean and the optimal estimate error of the i+1th iteration as the covariance matrix; D KL represents the KL divergence between two Gaussian distributions; γ represents a preset threshold; and i is an integer that increases from 0.

[0039] According to another aspect of the present disclosure, a Gaussian approximation filtering device for a nonlinear system of an industrial object is provided, comprising: a processor; a memory for storing instructions executable by the processor; wherein the processor is configured to implement the above method when executing the instructions stored in the memory.

[0040] According to another aspect of the present disclosure, a non-volatile computer-readable storage medium is provided, on which computer program instructions are stored, wherein the computer program instructions implement the above method when executed by a processor.

[0041] According to another aspect of the present disclosure, a computer program product is provided, including a computer-readable code, or a non-volatile computer-readable storage medium carrying the computer-readable code. When the computer-readable code runs in a processor of an electronic device, the processor in the electronic device executes the above method.

[0042] By mathematically modeling the physical system to be state estimated, a nonlinear state space model corresponding to the physical system is obtained; based on the nonlinear state space model, the first Gaussian distribution of the system state at time t is matched with a given distribution through moment matching operation to determine the predicted estimate and predicted estimate error of the system state at time t; the given distribution is based on the optimal estimate of the system state at time t-1 as the mean and the optimal estimate error of the system state at time t-1 as the covariance matrix; the system observation value at time t is obtained; based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, the optimal estimate and optimal estimate error of the system state at time t are iteratively updated until the optimal estimate and optimal estimate error at time t are obtained when the iteration stop condition is met; wherein the second Gaussian distribution takes the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix; on the one hand, linearization error can be avoided by minimizing the update cost in the update step instead of linearizing the model to approximate the steady-state point. At the same time, using natural gradient descent, utilizing the Fisher information matrix to adjust the gradient direction, and considering the curvature of the parameter space, combining the reconstruction update step with the moment matching in the prediction step can provide a more accurate system state estimation for nonlinear systems.

[0043] Further features and aspects of the present disclosure will become apparent from the following detailed description of exemplary embodiments with reference to the attached drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0044] The accompanying drawings, which are incorporated in and constitute a part of the specification, illustrate exemplary embodiments, features, and aspects of the disclosure and, together with the description, serve to explain the principles of the disclosure.

[0045] Figure 1 A flowchart showing a Gaussian approximation filtering method for a nonlinear system of an industrial object according to an embodiment of the present disclosure is shown;

[0046] Figure 2 A schematic diagram showing performance comparison results of different system state estimation methods according to an embodiment of the present disclosure;

[0047] Figure 3 A block diagram of a Gaussian approximation filtering device for a nonlinear system of an industrial object according to an embodiment of the present disclosure is shown;

[0048] Figure 4 A block diagram of a Gaussian approximation filtering device for a nonlinear system of an industrial object according to an embodiment of the present disclosure is shown. DETAILED DESCRIPTION

[0049] Various exemplary embodiments, features and aspects of the present disclosure will be described in detail below with reference to the accompanying drawings. The same reference numerals in the accompanying drawings represent elements with the same or similar functions. Although various aspects of the embodiments are shown in the accompanying drawings, the drawings are not necessarily drawn to scale unless otherwise specified.

[0050] The word “exemplary” is used exclusively herein to mean “serving as an example, example, or illustration.” Any embodiment described herein as “exemplary” is not necessarily to be construed as preferred or advantageous over other embodiments.

[0051] In addition, in order to better illustrate the present disclosure, numerous specific details are given in the following specific embodiments. It should be understood by those skilled in the art that the present disclosure can also be implemented without certain specific details. In some examples, methods, means, components and circuits well known to those skilled in the art are not described in detail in order to highlight the subject matter of the present disclosure.

[0052] In order to solve the technical problem that in the traditional system state estimation method, due to the approximate linearization of the nonlinear system, the estimation result introduces a large linearization error, the application reconstructs the update step of the Gaussian filter and proposes a new iterative update method, which directly minimizes the target of the update step, rather than linearizing the model to approximate the steady-state point, thereby avoiding the linearization error. At the same time, in the reconstructed update step, natural gradient descent is used, the Fisher information matrix is ​​used to adjust the gradient direction, and the curvature of the parameter space is considered. The reconstructed update step is combined with the moment matching in the prediction step, and a natural gradient Gaussian approximation filter (abbreviated as NANO filter) is introduced. The NANO filter can provide a more accurate system state estimation of the nonlinear system.

[0053] Among them, natural gradient descent is an optimization algorithm that improves the traditional gradient descent method by considering the geometric structure of the parameter space. In natural gradient descent, the direction of parameter update depends not only on the gradient of the loss function, but also on the correlation between parameters, which is usually achieved by using the Fisher information matrix.

[0054] The following is a detailed introduction to the Gaussian approximate filtering method for industrial object-oriented nonlinear systems provided by the present application.

[0055] Figure 1 The flowchart of the Gaussian approximate filtering method for nonlinear systems of industrial objects according to an embodiment of the present disclosure is shown. This embodiment takes the method used in an electronic device with computing capabilities as an example for explanation. The electronic device may be a user terminal or a server. The user terminal includes but is not limited to: a computer, a mobile phone, a tablet computer, etc. This embodiment does not limit the implementation method of the electronic device. Figure 1 As shown, the method includes:

[0056] Step 101 , mathematically modeling a physical system for which state estimation is to be performed to obtain a nonlinear state space model corresponding to the physical system.

[0057] Among them, the physical system to be state estimated refers to a dynamic system in reality, such as a vehicle control system, a power system, etc. The system state of the physical system changes with time, and this embodiment does not limit the implementation method of the physical system. In this embodiment, the physical system can be an industrial object with complex variables, nonlinearity, and / or variable environmental interference. In this embodiment, the system state refers to the system variables in the physical system, such as the system state is the position of the vehicle and / or the speed of the vehicle. Estimating the system state can realize monitoring and early warning of the physical system, etc. This embodiment does not limit the application scenario of the estimated results of the system state.

[0058] In one example, a physical system to be estimated is mathematically modeled, and the resulting nonlinear state space model can be expressed as follows:

[0059] x t+1 =f(x t )+ξ t ;

[0060]

[0061] in, represents the system state of the physical system, or in other words, x t is the real state or actual value of the physical system at time t. In this practical application, x t is an unknown variable, so we need to t The optimal estimate of the value of . represents n-dimensional Euclidean space. is the system observation value of the physical system, the system observation value y t It refers to the value that can be obtained by observing or measuring the physical system at time t. It is the data that can be directly obtained from the physical system. This value usually contains some noise The system observation value is the external manifestation of the system state, which is not directly equivalent to the system state x t , but the system observations can be used to update the system state x t Estimates. Represents m-dimensional Euclidean space, n and m are positive integers, and the values ​​of n and m are determined based on the estimation requirements of the system state. For example, if there are 5 system states to be estimated for the physical system, then the value of n is 5. represents the transfer matrix, and f is determined based on the dynamic model of the physical system. Represents the observation matrix, g is based on the system observations and the system state x to be estimatedt The conversion relationship between them is determined.

[0062] ξ t represents the process noise of the physical system; represents the observation noise of the physical system, ξ t and They are all approximately Gaussian distributions with a mean of 0, that is:

[0063]

[0064] Among them, 0 n×1 Represents a matrix with n rows and 1 column of elements equal to 0, 0 m×1 represents a matrix with m rows and 1 column with elements equal to 0. is the covariance matrix of process noise, is the covariance matrix of the observation noise; where Q and R are known matrices, representing the uncertainties of process noise and observation noise, respectively.

[0065] Process noise refers to the noise caused by random disturbances of the physical system, such as noise caused by model simplification, external interference of the physical system or other unknown factors. Q is a semi-positive definite matrix, that is, there may be some zero eigenvalues ​​in Q, which indicate the system state without random disturbances.

[0066] Observation noise refers to the random error in the observation process. Observation noise always exists and is unpredictable. Therefore, R is a positive definite matrix to indicate that all system observations have a certain degree of uncertainty.

[0067] Q and R can be determined based on empirical values ​​or by using parameter identification technology. This embodiment does not limit the method for determining Q and R.

[0068] ξ at different times t Independent of each other, ξ at different times t Independent of each other, and ξ at any time t and Independent of each other and independent of the initial state x0 of the physical system; the initial state x0 follows a Gaussian distribution in, and are the preset initial values ​​of random variables and the initial values ​​of the uncertainty of the initial state respectively.

[0069] Step 102, based on the nonlinear state space model, through moment matching operation, the first Gaussian distribution of the system state at time t is matched with a given distribution to determine the predicted estimate and predicted estimate error of the system state at time t; wherein the given distribution takes the optimal estimate of the system state at time t-1 as the mean and the optimal estimate error of the system state at time t-1 as the covariance matrix. t is a positive integer.

[0070] In the variational problem of Bayesian filtering, the system state x t The prior distribution of can be viewed as maximizing the variational probability function q(x t ) about the joint distribution p(x t-1 |y 1:t-1 ) and the transition probability p(x t |x t-1 ), namely:

[0071]

[0072] Among them, p(x t-1 |y 1:t-1 ) indicates that under the condition of the system observation value from time 1 to time t-1, the system state x is obtained t-1 The conditional probability distribution of the system state x t The prior distribution of p(x t |x t-1 ) represents the system state x at time t-1 t-1 Transition to the system state x at time t t The conditional probability distribution of ; It means that under the given condition p(x t-1 |y 1:t-1 ) and p(x t |x t-1 ), the variational probability function q(x t )’s expected value; It means to find the variational probability function q(x t ).

[0073] System Status x t The posterior distribution of can be viewed as minimizing the variational probability function q(x t ) about the observation probability p(y t |x t ) and the solution of the sum of the relative entropy (Kullback-Leibler divergence, KL divergence) between the variational probability function and the prior distribution, that is:

[0074]

[0075] p(x t |y 1:t ) indicates that under the condition of the system observation value from time 1 to time t, the system state x is obtained t The conditional probability distribution of the system state x t The posterior distribution of p(y t |x t) represents the system state x t Get the system observation value y t The conditional probability distribution of ; Indicates that in the variational probability function q(x t ) condition, for the system state x t The expected value of the negative log-likelihood of ; It means to find the variational probability function q(x) that minimizes the expected value of the negative log-likelihood. t ). KL Represents the KL divergence calculation function, which is used to measure q(x t ) and p(x t |y 1:t-1 )The difference between these two distributions, p(x t |y 1:t-1 ) indicates that under the condition of the system observation value from time 1 to time t-1, the system state x is obtained t The conditional probability distribution of .

[0076] In this embodiment, the variational probability function That is, the variational probability function q maps the n-dimensional real number space to a one-dimensional real number space.

[0077] According to the Gaussian approximation in the following formula:

[0078]

[0079] The above prior distribution can be described as:

[0080]

[0081] The posterior distribution is described as:

[0082]

[0083] in, represents the prediction cost, Represents the update cost. It means that the system state x at time t is known when the system observation value from time 1 to time t-1 is known. t The predicted estimate of the system state x t A priori estimate of P t|t-1 Represents the system state x t The forecast estimate error is used to indicate the uncertainty of the forecast estimate. Indicates that the system state x is known when the system observation value from time 1 to time t is known. t The posterior estimate of P t|t Represents the system state x tThe a posteriori estimation error is used to indicate the uncertainty of the a posteriori estimation. Represents the system state x t The best estimate of P t Represents the system state x t The optimal estimation error is used to indicate the uncertainty of the optimal estimation. It means that x follows a Gaussian distribution with μ as the mean and ∑ as the covariance matrix. The same symbol parameters have the same meaning as here and will not be repeated below.

[0084] The stable point of can be simplified to the moment matching equation. Specifically, for the random probability density function p(x), the stable point of the problem of maximizing the expected Gaussian likelihood can be expressed as follows:

[0085]

[0086] The above formula can be converted to:

[0087]

[0088] Among them, μ * To maximize the mean of the expected Gaussian likelihood stable point, ∑ * is the variance of the stable point of maximizing the expected Gaussian likelihood. According to the stable point model of the problem of maximizing the expected Gaussian likelihood, the stable point of maximizing the expected Gaussian likelihood is achieved when the Gaussian distribution of the variable x matches the mean and variance of the given distribution. Based on this, the system state x can be matched by the moment matching operation. t The first Gaussian distribution matches the given distribution Thus, the prior distribution obtained by the above Gaussian approximation is converted to:

[0089]

[0090] Among them, moment matching is a method used in statistics to estimate the parameters of probability distribution. It determines the value of the parameter by matching the moments of the sample (such as mean, variance, etc.) with the moments of the theoretical distribution. This method assumes that the sample data can represent the characteristics of the population, and estimates the parameters by making the sample moments as close as possible to the population moments. Moment matching can be used to estimate the parameters of various probability distributions, including Gaussian distribution, Poisson distribution, etc. In this embodiment, moment matching uses the given distribution as the moment of the sample and the system state x t As the parameters of the Gaussian distribution to be estimated, for the system state x t Make an estimate.

[0091] Based on the nonlinear state space model, the above prior distribution can be further transformed into:

[0092]

[0093] in, represents the predicted estimate of the system state at time t; P t|t-1 represents the prediction estimation error of the system state at time t; represents the given distribution; f represents the transfer matrix in the nonlinear state space model; ξ t-1 The matrix representing the process noise at time t-1 in the nonlinear state space model; represents the expected value calculation function; Represents the variance calculation function.

[0094] Therefore, the optimal estimate of the system state at time t-1 can be obtained by the above formula and the optimal estimation error P t-1|t-1 , determine the predicted estimate of the system state at time t and the prediction estimation error P t|t-1 .

[0095] The optimal estimate of the system state at initialization time 0 and the optimal estimation error P 0|0 is the preset value, that is, and All are preset known quantities.

[0096] Step 103, obtaining the system observation value at time t.

[0097] Schematically, data is collected or read from the physical system to obtain the system observation value at time t.

[0098] Step 104, based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, iteratively update the optimal estimate and the optimal estimate error of the system state at time t until the iteration stop condition is met, and the optimal estimate and the optimal estimate error at time t are obtained; wherein the second Gaussian distribution takes the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix.

[0099] In one example, based on a natural gradient descent method and a Fisher information matrix associated with a second Gaussian distribution of the system state at time t, iteratively updating an optimal estimate and an optimal estimate error of the system state at time t includes:

[0100] An iterative update model is provided for obtaining the optimal estimate and the optimal estimate error of the system state at time t; the predicted estimate and the predicted estimate error of the system state at time t, the optimal estimate and the optimal estimate error of the i-th iteration at time t, and the system observation value at time t are input into the iterative update model to obtain the optimal estimate and the optimal estimate error of the i+1-th iteration.

[0101] The iterative update model is obtained by stacking the Gaussian parameters of the second Gaussian distribution into a single column vector, determining the derivative with respect to the column vector, and then combining the column vector, the derivative with respect to the column vector, the Fisher information matrix, and the expression of the natural gradient descent method.

[0102] The following is an introduction to the acquisition process of the iterative update model. In order to simplify the iterative update model, in this embodiment, the Gaussian parameters of the second Gaussian distribution (i.e., the mean and covariance matrix of the second Gaussian distribution, where the mean is the optimal estimate of the system state at time t, and the covariance matrix is ​​the optimal estimate error of the system state at time t) are stacked into a single column vector v;

[0103]

[0104] Vec represents the matrix P t -1 Convert to a column vector, Represents the system state x t The best estimate of P t Represents the system state x t The optimal estimation error.

[0105] Calculating the derivative of the update cost with respect to the column vector yields:

[0106]

[0107] In this embodiment, the value of the variable v is adjusted by the update amount δv, and δv is expressed by the following formula:

[0108]

[0109] Here, i represents the number of iterations (or iteration index), and i is an integer that increases from 0.

[0110] Minimize the update cost by natural gradient descent The obtained update amount δv can be expressed by the following formula:

[0111]

[0112] in, Represents the second Gaussian distribution The inverse of the associated Fisher information matrix is ​​expressed in order to optimize the update cost Calculate the update cost using the value of v at iteration i The gradient of the variable v and the inverse of the Fisher information matrix give the update δv.

[0113] With the second Gaussian distribution The inverse of the associated Fisher information matrix is ​​expressed as follows:

[0114]

[0115] in, represents the Kronecker product.

[0116] Combining the above column vector, the update amount δv represented by the natural gradient descent method and the Fisher information matrix, the combined model is expressed as follows:

[0117]

[0118] in, represents the system state x obtained at the i+1th iteration at time t t The optimal estimation error The inverse of represents the system state x obtained at the i+1th iteration at time t t The best estimate of The system state x obtained in the i-th iteration t The optimal estimation error P t The inverse of represents the system state x obtained at the i-th iteration t The best estimate of .

[0119] According to Gaussian distribution Expectations The gradient with respect to the mean μ satisfies the following lemma:

[0120]

[0121] Gaussian distribution Expectations The Hessian matrix (i.e., the square matrix of second-order partial derivatives) of the mean μ satisfies the following lemma:

[0122]

[0123] In addition, the update cost On the optimal estimate of the system state and the optimal estimation error P tThe partial derivative of can be expressed as follows through the pre-built loss function:

[0124]

[0125] Combining the above lemma and the update cost About the Optimal Estimation and the optimal estimation error P t The partial derivative formula of can convert the combined model into a derivative-free form:

[0126]

[0127]

[0128] According to the above derivative-free formula, the iterative update model can be expressed as follows:

[0129]

[0130] in, represents the optimal estimation error of the i+1th iteration The inverse of represents the optimal estimate of the i+1th iteration; represents the prediction estimation error of the system state at time t; represents the predicted estimate of the system state at time t; represents the optimal estimation error of the i-th iteration at time t The inverse of represents the optimal estimate of the i-th iteration at time t; represents the preset loss function; represents the system state x at time t t It follows a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; y t represents the system observation value at time t; Represents the expected value calculation function. The same symbol parameters have the same meaning as here and will not be repeated in the following.

[0131] The optimal estimate and the optimal estimate error at the 0th iteration are set based on a preset method. In one example, in order to ensure that the covariance matrix in the iterative update model is positive definite, the preset method includes setting the optimal estimate and the optimal estimate error at the 0th iteration by Laplace approximation, and the preset method is expressed by the following formula:

[0132]

[0133] in, represents the optimal estimate at the 0th iteration at time t, represents the inverse of the optimal estimation error at the 0th iteration at time t; Represents the system state x estimated based on the maximum a posteriori probability at time t t The best estimate of represents the system state x at time t t It obeys the Gaussian distribution with the prediction estimate at time t as the mean and the prediction estimate error at time t as the covariance matrix; l represents the preset loss function; y t represents the system observation value at time t; exp represents the exponential operation function; Indicates the search for The system state x with the largest value t ; Indicates the system status Department, About system status x t The second derivative of .

[0134] In other embodiments, the preset method may also be Cholesky decomposition or square root parameterization to ensure the positive definiteness of the covariance matrix. This embodiment does not limit the implementation method of the preset method.

[0135] Optionally, the loss function includes but is not limited to one of the following:

[0136] 1. The observation likelihood function is expressed as follows:

[0137] l(x t ,y t )=-log p(y t |x t );

[0138] l(x t ,y t ) represents the loss function value output by the observation likelihood function, x t represents the system state at time t, y t represents the system observation value at time t. According to the expression of the observation likelihood function, the loss function represents the system state x t The system observation value y is observed under the condition t The smaller the negative log likelihood, the better the observation model p(y t |x t ) is more consistent with the actual observation results.

[0139] 2. The likelihood function induced by Beta divergence is expressed as follows:

[0140]

[0141] Among them, l β(x t ,y t ) represents the loss function value of the likelihood function output induced by Beta divergence, x t represents the system state at time t, y t represents the system observation value at time t. β is a preset value, for example, β = 0.01. The value of β can also be other values. This embodiment does not limit the value of β. According to the likelihood function induced by Beta divergence, it can be known that the loss function combines the observation model p(y t |x t ) can measure the accuracy of the observation model and control the sensitivity of the loss function to different terms through β.

[0142] 3. Pseudo-Huber loss function, expressed by the following formula:

[0143]

[0144] in, ph (x t ,y t ) represents the loss function value output by the pseudo-Huber loss function, δ represents the threshold variable, δ>0; x t represents the system state at time t, y t represents the system observation value at time t, and g represents the observation matrix in the nonlinear state space model.

[0145] 4. Weighted Observation Likelihood Function:

[0146] w (x t ,y t )=-w(x t ,y t )·logp(y t |x t );

[0147]

[0148] Among them, (x t ,y t ) represents the loss function value output by the weighted observation likelihood function, x t represents the system state at time t, y t represents the system observation value at time t, w: is the weight function, representing the n-dimensional real number space The vector sum comes from the m-dimensional real space The vector is mapped to the non-negative real space. c is a constant greater than 0; g represents the observation matrix in the nonlinear state space model; R is the observation noise in the nonlinear state space model The variance of .

[0149] In one example, the iteration stop condition is determined based on the divergence of the Gaussian distribution corresponding to two consecutive iterations. Schematically, the iteration stop condition is expressed by the following formula:

[0150]

[0151] in, represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i+1th iteration as the mean and the optimal estimate error of the i+1th iteration as the covariance matrix; D KL represents the KL divergence between two Gaussian distributions; γ represents the preset threshold; i is an integer that increases from 0.

[0152] When the Gaussian distribution corresponding to the i-th iteration and the Gaussian distribution corresponding to the i+1-th iteration meet the above iteration stopping condition, the iteration stops. At this time, the optimal estimate at time t is the optimal estimate obtained by the i-th iteration, that is, Correspondingly, the optimal estimation error at time t is the optimal estimation error obtained in the i-th iteration, that is,

[0153] In summary, the Gaussian approximate filtering method for nonlinear systems for industrial objects provided in this embodiment mathematically models the physical system to be state estimated to obtain a nonlinear state space model corresponding to the physical system; based on the nonlinear state space model, the first Gaussian distribution of the system state at time t is matched with the given distribution through moment matching operation to determine the predicted estimate and predicted estimate error of the system state at time t; the given distribution is based on the optimal estimate of the system state at time t-1 as the mean and the optimal estimate error of the system state at time t-1 as the covariance matrix; the system observation value at time t is obtained; based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, the optimal estimate and optimal estimate error of the system state at time t are iteratively updated until the optimal estimate and optimal estimate error at time t are obtained when the iteration stop condition is met; wherein the second Gaussian distribution takes the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix; on the one hand, the linearization error can be avoided by minimizing the update cost in the update step instead of linearizing the model to approximate the steady-state point. At the same time, using natural gradient descent, utilizing the Fisher information matrix to adjust the gradient direction, and considering the curvature of the parameter space, combining the reconstruction update step with the moment matching in the prediction step can provide a more accurate system state estimation for nonlinear systems.

[0154] In order to illustrate the robustness and accuracy of the Gaussian approximate filtering method for industrial nonlinear systems provided by this application, a simulation example is given below to illustrate the method. In this example, the physical model is an air traffic control scene. The aircraft performs a maneuver turn on a horizontal plane at an altitude of h = 50 relative to the radar. The state vector is defined as x = [p x ,v x ,p y ,v y ,ω] T Take 2 as an example to illustrate. x and p y represents the position of the aircraft, v x ,v y represents the speed of the aircraft in the plane coordinate system, and ω represents the turning rate of the aircraft.

[0155] The physical model can be converted into a nonlinear state space model as follows:

[0156]

[0157] According to the above embodiment, x t is the system state of the physical system to be estimated at time t, y t Refers to the system observation value at time t. ξ trepresents the process noise of the physical system, represents the observation noise of the physical system. In this example, Among them, Q and R are both known matrices, which are expressed by the following formulas:

[0158]

[0159] Among them, q1, q2, σ r , σ θ , Determined based on empirical values, for example: all values ​​recorded in literature, such as: q1 = 0.5, q2 = 1 × 10 -6 , Δt represents the sampling period, and Δt is a fixed value, for example: Δt=0.2.

[0160] In this example, the traditional extended Kalman filter (Extended Kalman Filter, EKF), unscented Kalman filter (Unscented Kalman Filter, UKF), iterative extended Kalman filter (IteratedExtended Kalman Filter, IEKF) and posterior linearization filter (Posterior LinearizationFilter, PLF) are compared with the state estimation method provided in this application. Through 100 Monte Carlo experiments, each experiment includes 200 time steps (that is, the maximum value of t in the above text is 200), the root mean square error (Root Mean Square Error, RMSE) of the optimal estimate of each time step is evaluated. The evaluation results are referenced Figure 2 As shown. Figure 2 It can be seen that in addition to estimating the RMSE of the optimal estimates obtained by EKF, UKF, IEKF and PLF, this example also estimates the RMSE of the optimal estimates obtained by the state estimation method provided by this application when using different loss functions and when each loss function uses different parameters. Specifically, the method using the observation likelihood function corresponds to Figure 2 The NANO legend in ; the method using the likelihood function induced by Beta divergence corresponds to Figure 2 In the β-NANO legend, the values ​​of β in the likelihood function induced by Beta divergence include 5×10 -3 , 1×10 -2 , 2×10 -2 Three methods; using pseudo-Huber loss function corresponds to Figure 2 In the Huber-NANO legend, the values ​​of δ in the pseudo-Huber loss function include 3, 5, and 10; the method using the weighted observation likelihood function corresponds to Figure 2In the Weight-NANO legend, the values ​​of c in the weighted observation likelihood function include 5, 10, and 15. Figure 2 It can be seen that the state estimation method provided in this application uses no matter which loss function (such as Figure 2 For any parameter in the state estimation method, the RMSE is smaller than the RMSE of the traditional state estimation method, that is, the state estimation method of the present application has better accuracy and anti-interference ability.

[0161] Figure 3 A block diagram of a Gaussian approximation filtering device for a nonlinear system of an industrial object according to an embodiment of the present disclosure is shown as follows: Figure 3 As shown, the device includes at least the following modules: a system modeling module 310 , a state prediction module 320 , a state observation module 330 and a state update module 340 .

[0162] A system modeling module 310 is used to perform mathematical modeling on a physical system for which state estimation is to be performed, and obtain a nonlinear state space model corresponding to the physical system;

[0163] The state prediction module 320 is used to match the first Gaussian distribution of the system state at time t with a given distribution through moment matching operation based on the nonlinear state space model to determine the predicted estimate and predicted estimate error of the system state at time t; the given distribution takes the optimal estimate of the system state at time t-1 as the mean and the optimal estimate error of the system state at time t-1 as the covariance matrix; the optimal estimate and the optimal estimate error of the system state at time 0 are initialized as preset values; t is a positive integer;

[0164] A state observation module 330, used to obtain the system observation value at time t;

[0165] The state update module 340 is used to iteratively update the optimal estimate and optimal estimate error of the system state at time t based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, until the optimal estimate and optimal estimate error at time t are obtained when the iteration stop condition is met; wherein the second Gaussian distribution has the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix.

[0166] For relevant details, please refer to the above method embodiment.

[0167] In some embodiments, the functions or modules included in the device provided by the embodiments of the present disclosure can be used to execute the method described in the above method embodiments. The specific implementation can refer to the description of the above method embodiments, and for the sake of brevity, it will not be repeated here.

[0168] The embodiment of the present disclosure also provides a computer-readable storage medium on which computer program instructions are stored, and the computer program instructions implement the above method when executed by a processor. The computer-readable storage medium can be a volatile or non-volatile computer-readable storage medium.

[0169] An embodiment of the present disclosure further proposes an electronic device, comprising: a processor; and a memory for storing instructions executable by the processor; wherein the processor is configured to implement the above method when executing the instructions stored in the memory.

[0170] The embodiments of the present disclosure also provide a computer program product, including a computer-readable code, or a non-volatile computer-readable storage medium carrying the computer-readable code. When the computer-readable code runs in a processor of an electronic device, the processor in the electronic device executes the above method.

[0171] Figure 4 1 is a block diagram of a Gaussian approximation filtering device 1900 for a nonlinear system of an industrial object according to an exemplary embodiment. For example, the device 1900 can be provided as a server or a terminal device. Figure 4 , the apparatus 1900 includes a processing component 1922, which further includes one or more processors, and a memory resource represented by a memory 1932 for storing instructions, such as an application, that can be executed by the processing component 1922. The application stored in the memory 1932 may include one or more modules, each corresponding to a set of instructions. In addition, the processing component 1922 is configured to execute instructions to perform the above method.

[0172] The device 1900 may also include a power supply component 1926 configured to perform power management of the device 1900, a wired or wireless network interface 1950 configured to connect the device 1900 to a network, and an input / output interface 1958 (I / O interface). The device 1900 may operate based on an operating system stored in the memory 1932, such as Windows Server 2000. TM , MacOS X TM , Unix TM ,Linux TM , FreeBSD TM or similar.

[0173] In an exemplary embodiment, a non-volatile computer-readable storage medium is also provided, such as a memory 1932 including computer program instructions, which can be executed by the processing component 1922 of the device 1900 to perform the above method.

[0174] The embodiments of the present disclosure have been described above, and the above description is exemplary, not exhaustive, and is not limited to the disclosed embodiments. Many modifications and changes will be apparent to those of ordinary skill in the art without departing from the scope and spirit of the described embodiments. The selection of terms used herein is intended to best explain the principles of the embodiments, practical applications, or technical improvements in the market, or to enable other persons of ordinary skill in the art to understand the embodiments disclosed herein.

Claims

1. A Gaussian approximate filtering method for nonlinear systems of industrial objects, characterized in that: The method comprises: Mathematically modeling the physical system for which state estimation is to be performed to obtain a nonlinear state space model corresponding to the physical system; Based on the nonlinear state space model, a first Gaussian distribution of the system state at time t is matched with a given distribution through a moment matching operation to determine a predicted estimate and a predicted estimate error of the system state at time t; the given distribution takes the optimal estimate of the system state at time t-1 as a mean and the optimal estimate error of the system state at time t-1 as a covariance matrix; the optimal estimate and the optimal estimate error of the system state at time 0 are initialized as preset values; and t is a positive integer; Get the system observation value at time t; Based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t, the optimal estimate and the optimal estimate error of the system state at time t are iteratively updated until the iteration stop condition is met, and the optimal estimate and the optimal estimate error at time t are obtained; wherein the second Gaussian distribution has the optimal estimate of the system state at time t as the mean and the optimal estimate error of the system state at time t as the covariance matrix.

2. The method according to claim 1, characterized in that: The method of iteratively updating the optimal estimate and the optimal estimate error of the system state at time t based on the natural gradient descent method and the Fisher information matrix associated with the second Gaussian distribution of the system state at time t includes: Obtaining an optimal estimate of the system state at the time t and an iterative update model of the optimal estimate error; wherein the iterative update model is obtained by stacking the Gaussian parameters of the second Gaussian distribution into a single column vector, determining the derivative with respect to the column vector, and then combining the column vector, the derivative with respect to the column vector, the Fisher information matrix, and the expression of the natural gradient descent method; The predicted estimate and predicted estimate error of the system state at time t, the optimal estimate and optimal estimate error of the i-th iteration at time t, and the system observation value at time t are input into the iterative update model to obtain the optimal estimate and optimal estimate error of the i+1-th iteration; wherein i is an integer increasing from 0, and the optimal estimate and optimal estimate error at the 0th iteration are set based on a preset method.

3. The method according to claim 2, characterized in that: The iterative update model is expressed by the following formula: in, represents the optimal estimation error of the i+1th iteration The inverse of represents the optimal estimate of the i+1th iteration; represents the prediction estimation error of the system state at time t; represents the predicted estimate of the system state at time t; represents the optimal estimation error of the i-th iteration at time t The inverse of represents the optimal estimate of the i-th iteration at time t; l represents the preset loss function; represents the system state x at time t t It follows a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; y t represents the system observation value at time t; Represents the expected value calculation function.

4. The method according to claim 2, characterized in that: The preset method includes setting the optimal estimate and the optimal estimate error at the 0th iteration by Laplace approximation, and the preset method is expressed by the following formula: in, represents the optimal estimate at the 0th iteration at time t, represents the inverse of the optimal estimation error at the 0th iteration at time t; Represents the system state x estimated based on the maximum a posteriori probability at time t t The best estimate of represents the system state x at time t t It obeys the Gaussian distribution with the prediction estimate at time t as the mean and the prediction estimate error at time t as the covariance matrix; l represents the preset loss function; y t represents the system observation value at time t; exp represents the exponential operation function; Indicates the search for The system state x with the largest value t ; Indicates the system status Department, About system status x t The second derivative of .

5. The method according to claim 3 or 4, characterized in that: The loss function is one of the following: an observation likelihood function, a Beta divergence-induced likelihood function, a pseudo-Huber loss function, and a weighted observation likelihood function.

6. The method according to claim 1, characterized in that: The first Gaussian distribution of the system state at time t is matched to a given distribution through a moment matching operation based on the nonlinear state space model to determine the predicted estimate and the predicted estimate error of the system state at time t, which is expressed by the following formula: in, represents the predicted estimate of the system state at time t; P t|t-1 represents the prediction estimation error of the system state at time t; represents the given distribution; f represents the transfer matrix in the nonlinear state space model; ξ t-1 A matrix representing the process noise at time t-1 in the nonlinear state-space model; represents the expected value calculation function; Represents the variance calculation function.

7. The method according to claim 6, characterized in that: The nonlinear state space model is expressed by the following formula: x t+1 =f(x t )+ξ t ; y t =g(x t )+ζ t ; in, represents the system state of the physical system, represents n-dimensional Euclidean space; is the system observation value of the physical system, represents m-dimensional Euclidean space, n and m are positive integers; f: represents the transfer matrix; g: represents the observation matrix; ξ t represents the process noise of the physical system; t represents the observation noise of the physical system, ξ t and t They are all approximately Gaussian distributions with a mean of 0, that is: Among them, 0 n×1 Represents a matrix with n rows and 1 column of elements equal to 0, 0 m×1 represents a matrix with m rows and 1 column with elements equal to 0. is the covariance matrix of process noise, is the covariance matrix of the observation noise; ξ at different times t Independent of each other, ξ at different times t Independent of each other, and ξ at any time t and t are independent of each other and of the initial state x0 of the physical system; the initial state x0 follows a Gaussian distribution in, and are the preset initial values ​​of random variables and the initial values ​​of the uncertainty of the initial state respectively.

8. The method according to claim 1, characterized in that: The iteration stop condition is expressed by the following formula: in, represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i-th iteration as the mean and the optimal estimate error of the i-th iteration as the covariance matrix; represents the system state x at time t t It obeys a Gaussian distribution with the optimal estimate of the i+1th iteration as the mean and the optimal estimate error of the i+1th iteration as the covariance matrix; D KL represents the KL divergence between two Gaussian distributions; γ represents a preset threshold; and i is an integer that increases from 0.

9. A Gaussian approximate filtering device for nonlinear systems of industrial objects, characterized in that: include: processor; a memory for storing processor-executable instructions; The processor is configured to implement the method described in any one of claims 1 to 8 when executing the instructions stored in the memory.

10. A non-volatile computer-readable storage medium having computer program instructions stored thereon, characterized in that: When the computer program instructions are executed by a processor, the method according to any one of claims 1 to 8 is implemented.

Images