A Measurement-Data-Driven Control Approach towards Variance Reduction of Micromachined Resonant Accelerometer under Multi Unknown Disturbances

This paper first presents an adaptive expectation-maximization (AEM) control algorithm based on a measurement-data-driven model to reduce the variance of microelectromechanical system (MEMS) accelerometer sensor under multi disturbances. Significantly different characteristics of the disturbances, consisting of drastic-magnitude, short-duration vibration in the external environment, and slowly-varying, long-duration fluctuation inside the sensor are first constructed together with the measurement model of the accelerometer. Next, through establishing a data-driven model based on a historical small measurement sample, the window length of filter of the presented algorithm is adaptively chosen to estimate the sensor state and identify these disturbances simultaneously. Simulation results of the proposed AEM algorithm based on experimental test are compared with the Kalman filter (KF), least mean square (LMS), and regular EM (REM) methods. Variances of the estimated equivalent input under static condition are 0.212 mV, 0.149 mV, 0.015 mV, and 0.004 mV by the KF, LMS, REM, and AEM, respectively. Under dynamic conditions, the corresponding variances are 35.5 mV, 2.07 mV, 2.0 mV, and 1.45 mV, respectively. The variances under static condition based on the proposed method are reduced to 1.9%, 2.8%, and 27.3%, compared with the KF, LMS, and REM methods, respectively. The corresponding variances under dynamic condition are reduced to 4.1%, 70.1%, and 72.5%, respectively. The effectiveness of the proposed method is verified to reduce the variance of the MEMS resonant accelerometer sensor.


Introduction
Microelectromechanical system (MEMS) technologies have successfully enabled the miniaturization and cost reduction of many kinds of sensors. These MEMS sensors have been widely used in many fields [1][2][3], including artificial intelligence (AI), internet of things (IoT), and industry 4.0. Especially, the MEMS resonant accelerometer, as a kind of functional but important device, has recently focused lots of research interest due to its potentially high accuracy and built-in digital characteristic [4,5]. However, currently, the comparatively low performance of this type of accelerometer severely limits the application range, which suffers from lots of factors, such as scale-factor errors, bias drift, and additively harmful pressure and shock effects [6,7]. Among these factors, the comparatively large bias drift of the resonant accelerometer is one major reason to deteriorate its variance. Especially during a sharp vibration and shock occurrence, the bias of the accelerometer sensor will fluctuate dramatically, which leads to a terrible performance deterioration due to its sophisticated, but fragile mechanical structure [8]. Even worse, it is not easy to suppress these disturbances by common compensation methods, such as mechanical topology design [9][10][11], closed-loop circuit [5,12,13], and package technology [14][15][16] due to unknown and uncertain characteristic of these harmful inputs.
Signal post-processing technology is an effective way to estimate the unknown inputs in some fields, such as sensor modelling [17][18][19], object tracking [20][21][22], and coupled audio-visual analysis [23]. The bias model of MEMS sensors with unknown inputs (UIs) caused by inaccurate coefficient values, parameter variations, and uncertain disturbances can be estimated by establishing adaptive or robust controllers based on this post-processing technology. For example, a reduced-order disturbance decoupled observer is designed based on a parametric approach to jointly estimate the sensing system state and input with UI [24,25]. Further, a parameter-dependent robust mixed full-and reduced-order controllers for discrete-time sensing system with UI are obtained via the linear matrix inequality optimization [26]. UI of sensors are also modelled as a randomly switching parameter obeying a known Markov chain in multiple model estimation [27], where the corresponding controller contains variable-structure MME, interacting MME, and pseudo-Bayesian estimation. Another approach based on minimum upper bound filter (MUBF) is used to estimate a statistically-constrained UI by establishing the minimized the upper bounds of covariance matrices of the state prediction error and residual error [28].
For the MEMS resonant accelerometer sensor, the corresponding variance is usually affected by two types of inputs. The one, is vibration and shock with a characteristic of drastic amplitude and short duration sourcing from external environment. The other one, is a slowly varying fluctuation sourcing from the accelerometer's interior, such as thermal, pressure drift, and stress release of mechanical components. These two inputs can be commonly regarded as unknown parameters with stochastic and uncertain characteristics.
Currently, the data-driven-based signal processing approach is becoming a potentially interesting research to estimate the uncertain input of sensors [29][30][31][32][33]. The statistical models of sensors are fitted upon the available historical and condition data. In this work, the adaptive expectation-maximization (EM) control algorithm based on the data-driven model is first proposed to reduce the variance of the MEMS resonant accelerometer sensor under multi unknown disturbances. Two types of the external and internal disturbances are first considered simultaneously. Through constructing a data-driven model based on historical small sample data, a smoother window-length is adaptively chosen. The external and internal disturbances with significantly different characteristics are identified and the bias of the accelerometer is finally estimated. This paper is organized as follows: In Section 2, a model of the MEMS resonant accelerometer sensor under multi unknown disturbances is established. In Section 3, the joint estimation and identification algorithm based on the measurement-data-driven AEM are estimated and simulated. In Section 4, the experimental test of the proposed AEM method is compared and verified. Finally, the conclusions are given in Section 5.

Modelling of the Accelerometer Output Under Unknown Input
A typical MEMS resonant accelerometer sensor consists of a mechanical element and electrical element, as shown in Figure 1. y represents the output of MEMS accelerometer. When unknown inputs including external disturbance and internal disturbance occur, discrete-time model of output y k at time point t k can be first expressed as follows: where s k represents the true equivalent acceleration input of the MEMS accelerometer sensor. u ex,k represents the equivalent input induced by an external unknown disturbance at time point t k with the characteristic of drastic amplitude and short duration, such as shock and vibration. u in,k represents the equivalent input caused by an internal unknown disturbance at time point t k , such as stress release of the mechanical element and the thermal drift of electrical element components. w k represents zero-mean white Gaussian noise. Here, let u k = u ex,k + u in,k . where sk represents the true equivalent acceleration input of the MEMS accelerometer sensor. uex,k represents the equivalent input induced by an external unknown disturbance at time point tk with the characteristic of drastic amplitude and short duration, such as shock and vibration. uin,k represents the equivalent input caused by an internal unknown disturbance at time point tk, such as stress release of the mechanical element and the thermal drift of electrical element components. wk represents zero-mean white Gaussian noise. Here, let uk = uex,k + uin,k . A state-space model for the MEMS accelerometer sensor is first constructed. In most applications of MEMS sensors, the maneuverability characteristics of aircrafts can be regarded as constraining in a certain frequency bandwidth and magnitude. Additionally, the true input signal sk+1 of the sensors for the subsequent time point tk+1 usually can be considered as relating to sk for the former time point tk. The above two characteristics can be suitably described by a first-order Markov process with a process time constant set according to the system's bandwidth and the process noise related to the limit of the magnitude of the motion [34][35][36]. Hence, according to Equation (1), the dynamical and measurement models of the true input signal for MEMS accelerometer can be constructed, respectively. is described by a covariance matrix Q = σQ 2 I1×1 that can be set according to the dynamic characteristic of the true input signal. σQ 2 and σR 2 are the corresponding variance parameters, respectively. Coefficient scalar Φk = −1/τs, where τs is process time constant based on first-order Markov process. Hk+1, Γk+1 and ψk+1 are all one-dimensional identity matrixes. Here, symbol Θ represents vector {uex, uin, σQ 2 , σR 2 }.

Joint Estimation and Identification Strategy Based on the Data-Driven Iterative Optimization
To perform the identification of the unknown input uex, uin and the estimation of the true acceleration input signal s simultaneously, it is necessary to use an iterative optimization strategy to achieve joint estimation and identification because unknown input identification mistake deteriorates the state estimate while the state estimate error increases the identification risk. The EM algorithm is an effective approach to implement this iterative optimization. Nevertheless, due to the significant differences in the characteristics of magnitude and duration between the external and internal disturbances, it is not suitable for the EM algorithm to keep the window length of its filter at A state-space model for the MEMS accelerometer sensor is first constructed. In most applications of MEMS sensors, the maneuverability characteristics of aircrafts can be regarded as constraining in a certain frequency bandwidth and magnitude. Additionally, the true input signal s k+1 of the sensors for the subsequent time point t k+1 usually can be considered as relating to s k for the former time point t k . The above two characteristics can be suitably described by a first-order Markov process with a process time constant set according to the system's bandwidth and the process noise related to the limit of the magnitude of the motion [34][35][36]. Hence, according to Equation (1), the dynamical and measurement models of the true input signal for MEMS accelerometer can be constructed, respectively.
where k is the sampling time, x k = s k represents the state of accelerometer, z k+1 represents the measurement of sensor, u ex,k+1 and u in,k+1 are unknown external and internal disturbance, respectively. Similarly with Equation (1), is characterized by the white noise with covariance matrix R = σ R 2 I 1×1 that can be obtained by calculating the final output data of the single sensor based on the Allan variance method. δ(τ) is the Dirac delta function.
to the dynamic characteristic of the true input signal. σ Q 2 and σ R 2 are the corresponding variance parameters, respectively. Coefficient scalar Φ k = −1/τ s , where τ s is process time constant based on first-order Markov process. H k+1 , Γ k+1 and ψ k+1 are all one-dimensional identity matrixes. Here, symbol Θ represents vector {u ex , u in , σ Q 2 , σ R 2 }.

Joint Estimation and Identification Strategy Based on the Data-Driven Iterative Optimization
To perform the identification of the unknown input u ex , u in and the estimation of the true acceleration input signal s simultaneously, it is necessary to use an iterative optimization strategy to achieve joint estimation and identification because unknown input identification mistake deteriorates the state estimate while the state estimate error increases the identification risk. The EM algorithm is an effective approach to implement this iterative optimization. Nevertheless, due to the significant differences in the characteristics of magnitude and duration between the external and internal disturbances, it is not suitable for the EM algorithm to keep the window length of its filter at a constant. Hence, the AEM algorithm driven by measurement data of MEMS accelerometer sensor under complex disturbances is first presented to achieve iterative optimization, and the following derivation follows the idea of the proposed data-driven-based AEM strategy.
The function flow of the proposed method is shown in Figure 2. The proposed AEM method mainly consists of a likelihood function estimation, variable window-length smoother selection, state estimation in E step, parameter identification in M step and iteration decision. The detailed proposed method of the MEMS accelerometer sensor is given as follow.
The likelihood function of measurement z k of the MEMS accelerometer is first given to establish the corresponding maximum expectation. Then, in E-step and M-step stages, the maximum expectation of the log likelihood function is performed by these two iteration steps. The unknown external and internal disturbances are then identified based on their difference on window lengths of filters, which are driven by slope of acquired measurement data over a certain period of time. Finally, equivalent acceleration input of the sensor is estimated.
Micromachines 2019, 10, x FOR PEER REVIEW 4 of 12 a constant. Hence, the AEM algorithm driven by measurement data of MEMS accelerometer sensor under complex disturbances is first presented to achieve iterative optimization, and the following derivation follows the idea of the proposed data-driven-based AEM strategy. The function flow of the proposed method is shown in Figure 2. The proposed AEM method mainly consists of a likelihood function estimation, variable window-length smoother selection, state estimation in E step, parameter identification in M step and iteration decision. The detailed proposed method of the MEMS accelerometer sensor is given as follow.
The likelihood function of measurement zk of the MEMS accelerometer is first given to establish the corresponding maximum expectation. Then, in E-step and M-step stages, the maximum expectation of the log likelihood function is performed by these two iteration steps. The unknown external and internal disturbances are then identified based on their difference on window lengths of filters, which are driven by slope of acquired measurement data over a certain period of time. Finally, equivalent acceleration input of the sensor is estimated.

Likelihood Function of the Accelerometer Sensor Measurement
For the maximum-likelihood estimation of the MEMS accelerometer sensor, the complete-data log-likelihood function (CLLF) is defined to be the logarithm of the probability density function of output observations parameterized with Θ. First, the complete-data LLF for the parameter Θ can be written as follows: where v is the window length of the filter, p is the probability density function of the output measurement of the MEMS accelerometer. The function p of the above expression is assumed to obey Gaussian distribution and is expressed as follows:

Likelihood Function of the Accelerometer Sensor Measurement
For the maximum-likelihood estimation of the MEMS accelerometer sensor, the complete-data log-likelihood function (CLLF) is defined to be the logarithm of the probability density function of output observations parameterized with Θ. First, the complete-data LLF for the parameter Θ can be written as follows: where v is the window length of the filter, p is the probability density function of the output measurement of the MEMS accelerometer. The function p of the above expression is assumed to obey Gaussian distribution and is expressed as follows: Further, the log-likelihood function L k k−v+1 is given in detail as follows: with the following: where the dimensions d 1 and d 2 of the state x and measurement z are equal to 1, respectively. To estimate the state and identify the parameter of the MEMS accelerometer, maximum expectation of the log likelihood function is performed by the proposed AEM algorithm below.

Joint Estimation and Identification based on Iterative Optimization
The AEM algorithm to identify the disturbances and estimate the above equivalent input of MEMS accelerometer is first given as follows: where r is the number of iterations, The corresponding maximum expectation value is represented as Q r max because these partial derivations should be zero at the optimal point of state x estimation and disturbance u identification; therefore: Here, u i is equal to each other within the window length of v as follows: Therefore, we obtain the following: Then, the (r + 1) th -cycle Θ k|r+1 k−v+1 is iterated by replacing u k|r in Θ k|r k−v+1 with the above u k|r+1 i k−v+1 . Iteration will stop when any one of the below two cases occurs. The one is the number of iterations reaches the upper limit r max set by the simulation. The other one is the increment of the likelihood function between two consecutive iterations is less than the limit value δ, which means the following: where the range of δ is usually set as (0, 1), Q(r + 1) represents Q Θ k k−v+1 |Θ k|r+1 k−v+1 . To take the balance between output accuracy with the running speed of the algorithm, here δ is 0.01.

Adaptive Selection on Window Length of Filter
Due to the different characteristics of u ex and u in , the window-length of the filter is adaptively chosen. The larger smoothing window length of the filter is suitable to suppress the long-time slowly varying disturbance of the sensors effectively, such as u in , while the small window length is more effective to reduce short-time drastic fluctuation, such as u ex . In detail, slope variation criterion driven by historical small sample data is constructed to determine the window length of the filter in order to suppress these two kinds of disturbances simultaneously. Here, the window length model of the filter is formulated as follows: where ∆k s = max k s {z i ', z j '} − min k s {z i ', z j '} represents slope difference between the maximum slope value and the minimum value, k s {z i ', z j '} represents the slope value between any two points on the curve fitted based on historical sample measurement data. z i ' and z j ' are the new fitted sample points corresponding to historical sample points z i and z j , respectively. The recursive computation process of the derived AEM algorithm is summarized in Table 1. Table 1. Recursive computation process of the proposed method.
Step Description

2). Start
Set window length v of the filter

4). E-Step
Use KF to estimatex k|r k−v+1 , then calculate the expectation value E x k|r

5). M-Step
Identify the parameters u k|r 1 +1 k−v+1 to ensure the above expectation value to reach the maximum value Q r max

6). Termination
Calculate ∆Q in Equation (20), If ∆Q < δ or r > r max , then terminate the current number of iterations as r. Else, update Θ k|r k−v+1 as Θ k|r+1 k−v+1 by replacing u k|r k−v+1 with the above u k|r+1 k−v+1 , and go to step 4.

Experimental Test and Comparison
The test system of the MEMS sensor is mainly composed of a resonant accelerometer device [37], control circuit, data acquisition unit using Agilent 34410 with sampling frequency of 15 Hz to collect the output of sensors, and the PC that is connected to 34410 by USB-port communication interface and online process the measurement data by the proposed method. Figure 3 shows the overall test system of the MEMS sensor.

Static Performance Test
The MEMS resonant accelerometer sensor is performed to respond to the static output under stochastic disturbances with two types of characteristics of slowly varying fluctuation and drastic-amplitude and short-duration vibration. The results of the proposed method are compared with other vibration suppression algorithms, such as the LMS algorithm, the KF algorithm and the REM algorithm. The variance parameters are initially set as σQ = 0.93 and σR = 0.99.
First, some groups of stochastic signals, including a slowly varying fluctuation and drastic-amplitude, short-duration vibration are generated randomly. Figure 4 shows the estimated equivalent input by different algorithms. First, the detailed characteristics of the three groups of the disturbances are listed in Table 2. The equivalent weight of movable platform of the vibration table is 2 kg and the corresponding impact accelerometer is 490 m/s 2 . The corresponding results are listed in Table 3. The effectiveness of the proposed algorithm is proven by the variance of the equivalent input of 0.004 mV, compared with variances of 0.212 mV, 0.149 mV, and 0.015 mV estimated by KF, LMS, and REM algorithms, respectively. The corresponding variance by the proposed method is reduced to 1.9%, 2.8%, and 27.3%, respectively.

Static Performance Test
The MEMS resonant accelerometer sensor is performed to respond to the static output under stochastic disturbances with two types of characteristics of slowly varying fluctuation and drastic-amplitude and short-duration vibration. The results of the proposed method are compared with other vibration suppression algorithms, such as the LMS algorithm, the KF algorithm and the REM algorithm. The variance parameters are initially set as σ Q = 0.93 and σ R = 0.99.
First, some groups of stochastic signals, including a slowly varying fluctuation and drastic-amplitude, short-duration vibration are generated randomly. Figure 4 shows the estimated equivalent input by different algorithms. First, the detailed characteristics of the three groups of the disturbances are listed in Table 2. The equivalent weight of movable platform of the vibration table is 2 kg and the corresponding impact accelerometer is 490 m/s 2 .

Static Performance Test
The MEMS resonant accelerometer sensor is performed to respond to the static output under stochastic disturbances with two types of characteristics of slowly varying fluctuation and drastic-amplitude and short-duration vibration. The results of the proposed method are compared with other vibration suppression algorithms, such as the LMS algorithm, the KF algorithm and the REM algorithm. The variance parameters are initially set as σQ = 0.93 and σR = 0.99.
First, some groups of stochastic signals, including a slowly varying fluctuation and drastic-amplitude, short-duration vibration are generated randomly. Figure 4 shows the estimated equivalent input by different algorithms. First, the detailed characteristics of the three groups of the disturbances are listed in Table 2. The equivalent weight of movable platform of the vibration table is 2 kg and the corresponding impact accelerometer is 490 m/s 2 . The corresponding results are listed in Table 3. The effectiveness of the proposed algorithm is proven by the variance of the equivalent input of 0.004 mV, compared with variances of 0.212 mV, 0.149 mV, and 0.015 mV estimated by KF, LMS, and REM algorithms, respectively. The corresponding variance by the proposed method is reduced to 1.9%, 2.8%, and 27.3%, respectively. (a)

Dynamic Performance Test
The dynamic characteristic of the MEMS resonant accelerometer under the disturbances is also compared with other fluctuation filtering algorithms to verify the effectiveness of the proposed  The corresponding results are listed in Table 3. The effectiveness of the proposed algorithm is proven by the variance of the equivalent input of 0.004 mV, compared with variances of 0.212 mV, 0.149 mV, and 0.015 mV estimated by KF, LMS, and REM algorithms, respectively. The corresponding variance by the proposed method is reduced to 1.9%, 2.8%, and 27.3%, respectively.

Dynamic Performance Test
The dynamic characteristic of the MEMS resonant accelerometer under the disturbances is also compared with other fluctuation filtering algorithms to verify the effectiveness of the proposed method. The input acceleration has the amplitude and frequency of 2 mV and 0.2 Hz with a peak-to-peak value 0.11 mV of white Gaussian noise, respectively.
Further, the disturbances of the impulse and squared waves are applied to the resonant accelerometer randomly. The maximum amplitudes of these impulse signals are within the equivalent input voltage of ±60 mV. The detailed characteristics of the three groups of the squared wave and impulse disturbances are listed in Table 4. The equivalent input is estimated by different algorithms, as shown in Figure 5. method. The input acceleration has the amplitude and frequency of 2 mV and 0.2 Hz with a peak-to-peak value 0.11 mV of white Gaussian noise, respectively. Further, the disturbances of the impulse and squared waves are applied to the resonant accelerometer randomly. The maximum amplitudes of these impulse signals are within the equivalent input voltage of ±60 mV. The detailed characteristics of the three groups of the squared wave and impulse disturbances are listed in Table 4. The equivalent input is estimated by different algorithms, as shown in Figure 5.   The variance of the equivalent input estimation is listed in Table 5. The variance of the proposed algorithm is much smaller than the variances of other methods, which is decreased to 4.1%, 70.1%, and 72.5%, respectively. The effectiveness of the proposed method is proven.

Conclusions
The AEM control algorithm driven by measurement data is first established to reduce the variance of the MEMS resonant accelerometer sensor under multi unknown disturbances. The corresponding simulations based on experimental test are analyzed by comparing the proposed method with the KF, the LMS, and the REM methods. The variances of the equivalent input estimation of MEMS accelerometer with multi disturbances under static and dynamic conditions are 0.004 mV and 1.45 mV, respectively. The variances under static condition based on the proposed method are reduced to 1.9%, 2.8%, and 27.3%, compared with KF, LMS, and REM method, respectively. The corresponding variances under dynamic condition are reduced to 4.1%, 70.1%, and 72.5%, respectively. The effectiveness of the proposed method is verified to reduce the variance of the MEMS resonant accelerometer sensor.
Author Contributions: Q.S. conceived the idea of the paper, designed and performed the simulations and experiments, and wrote the manuscript; D.Y. did the simulations and analysis; J.Z. did the theoretical analysis; Y.W. did the review; Y.Z. did the visualization, and W.Y. reviewed the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.